Next Article in Journal
Assessment of Cavitation Models for Compressible Flows Inside a Nozzle
Next Article in Special Issue
Influence of Wind Buffers on the Aero-Thermal Performance of Skygardens
Previous Article in Journal / Special Issue
Numerical Simulation of Fountain Formation due to Normal and Inclined Twin-Jet Impingement on Ground
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

A Simple Construction of a Thermodynamically Consistent Mathematical Model for Non-Isothermal Flows of Dilute Compressible Polymeric Fluids

1
Faculty of Mathematics and Physics, Charles University, Sokolovská 83, CZ 186 75 Praha 8–Karlín, Czech Republic
2
Mathematical Institute, University of Oxford, Radcliffe Observatory Quarter, Woodstock Road, Oxford OX2 6GG, UK
*
Author to whom correspondence should be addressed.
Fluids 2020, 5(3), 133; https://doi.org/10.3390/fluids5030133
Submission received: 13 July 2020 / Revised: 4 August 2020 / Accepted: 6 August 2020 / Published: 11 August 2020
(This article belongs to the Special Issue Recent Advances in Fluid Mechanics: Feature Papers)

Abstract

:
We revisit some classical models for dilute polymeric fluids, and we show that thermodynamically consistent models for non-isothermal flows of these fluids can be derived in a very elementary manner. Our approach is based on the identification of energy storage mechanisms and entropy production mechanisms in the fluid of interest, which, in turn, leads to explicit formulae for the Cauchy stress tensor and for all of the fluxes involved. Having identified these mechanisms and derived the governing equations, we document the potential use of the thermodynamic basis of the model in a rudimentary stability analysis. In particular, we focus on finite amplitude (nonlinear) stability of a stationary spatially homogeneous state in a thermodynamically isolated system.
MSC:
76A05; 35Q79; 37L15

1. Introduction

Starting from the seminal work by Kramers [1] kinetic-type models have been widely used in the mathematical modelling of polymeric fluids, see the monographs by Bird et al. [2], Beris and Edwards [3], Öttinger [4], Huilgol and Phan-Thien [5], Dressler et al. [6], Öttinger [7] Kröger [8], and the review paper by Lozinski et al. [9] to name a few (a mathematically inclined reader—a specialist in theory of partial differential equations—is also referred to the concise presentation in Le Bris and Lelievre [10]). The overwhelming majority of the works, especially in the field of numerical simulations, see, for example Lozinski and Chauviere [11], Knezevic and Süli [12] or Mizerová and She [13], are however restricted to isothermal flows. In particular, the temperature field is often tacitly assumed to be homogeneous in space, and an evolution equation describing the temporal and spatial variations of the temperature field is rarely formulated, albeit some approaches, such as the GENERIC formalism, see Öttinger [7], or Pavelka et al. [14], allow one to do so. (regarding the classical macroscopic models, temperature evolution equations have been formulated for example in Peters and Baaijens [15], Wapperom and Hulsen [16] and Dressler et al. [6], see also Hron et al. [17] for further discussion). In the present contribution, we provide a straightforward self-contained derivation of a simple kinetic-type model for non-isothermal flows of compressible dilute polymeric fluids. The model is essentially the same as the model discussed in Öttinger and Grmela [18], but the approach to the derivation of a thermodynamically consistent model is substantially different.
We only use elementary arguments easily accessible to everyone familiar with basic principles of continuum mechanics such as the general form of balance equations. It turns out that, once one correctly specifies the energy storage mechanisms and the entropy production mechanisms in the fluid of interest, then some formulae, in particular the Kramers formula for the stress tensor, are natural consequences of the basic principles (there is no need to assume these formulae a priori). The same holds true in the derivation proposed by Öttinger and Grmela [18], where the authors claim that this fact is a “quite remarkable” consequence of the GENERIC structure. While the GENERIC framework has definitely its merits, we show that much more elementary arguments lead—regarding this issue—to the same conclusion. Furthermore, having obtained a thermodynamically consistent model, we exploit thermodynamic considerations in a rudimentary stability analysis in the spirit of Coleman [19], see also Bulíček et al. [20] for a discussion.
The paper is organised, as follows. First we recall basic principles of continuum mechanics, see Section 2, with particular attention placed on the basic tenets of the classical kinetic-type theory for dilute polymeric fluids. The key assumption regarding dilute polymeric fluids is that the polymeric chains do not effectively contribute to the mass of the solvent/polymer mixture, which means that the density of the mixture coincides with the density of the solvent. Furthermore, the kinetic energy of the polymeric chains is also neglected, the only mechanical energy contribution of the polymeric chains being via the energy stored in the stretched polymeric chains.
In Section 3, we focus on the dynamics of polymeric chains. Regarding the chains, we do not model them individually, but we again follow a kinetic-type theory, and we instead formulate a general form of the Fokker–Planck equation for the configurational distribution function of the polymeric chains. In our setting, the Fokker–Planck equation contains two unknown flux terms, namely the flux in configurational space and the flux in physical space. The latter is essential if we want to model the so-called centre-of-mass diffusion. Furthermore, we also comment on possible stationary solutions of the Fokker–Planck equation, and we carefully discuss the interplay between boundary conditions in the configurational space and asymptotic behaviour of the spring potential.
In Section 4, we give a formula for the specific Helmholtz free energy for the given fluid, which is we specify possible energy storage mechanisms. Subsequently, we proceed in Section 5 with the derivation of the constitutive relations, that is with the derivation of the equations relating the Cauchy stress tensor, the energy flux, the fluxes in the Fokker–Planck equation, and the kinematical quantities. The basic idea is that we look for closure relations that guarantee the nonnegativity of the entropy production. In this sense, the identification of the entropy production mechanisms and energy storage mechanisms provides a complete characterisation of the fluid of interest. A summary of the resulting model is given in Section 5.4.
Once we obtain the constitutive relations we discuss, see Section 6, the implications of the thermodynamic basis regarding the investigation of the stability of container flows (from the thermodynamic point of view, the flow in a thermally and mechanically isolated container constitutes an example of a thermodynamically isolated system). We show that the system of governing equations has a “trivial” stationary solution, and using thermodynamic arguments, we explicitly construct a nonnegative functional that decays in time and vanishes if and only if the system reaches the spatially homogeneous stationary state. The construction of such a functional is clearly a precursor for a rigorous stability analysis, which is, however, beyond the scope of the present contribution.

2. Preliminaries

The classical balance equations for the mass and momentum and the evolution equation for the internal energy of a continuous medium are
d ρ d t + ρ div x v = 0 ,
ρ d v d t = div x T + ρ b ,
ρ d e d t = T : D div x j e .
The symbol ρ denotes the density of the medium, v denotes the Eulerian velocity field, D = def 1 2 x v + x v denotes the symmetric velocity gradient, T denotes the Cauchy stress tensor, e denotes the specific internal energy, j e denotes the energy flux, and b denotes the external body force. The symbol
d d t = def t + v · x ,
denotes the material derivative, and the symbol A : B = def Tr A B denotes the Frobenius norm. The subscript x reminds us that the differential operator concerned is applied with respect to the spatial variable x . Regarding the derivation of the balance equations, see, for example, Málek and Průša [21] or any standard textbook on continuum mechanics, such as Müller [22] or Gurtin et al. [23].
The basic idea in a kinetic-type theory of dilute polymeric fluids is that the mass contribution of the polymer chains is assumed to be negligible, hence we can replace the total density ρ in (1) by the solvent density ρ s . Moreover, the barycentric velocity of the solvent/polymer mixture coincides with the velocity of the solvent, hence v in (1) can be interpreted as the velocity of the solvent only. The only place where the polymeric chains enter the balance equations is the formula for the internal energy e, the Cauchy stress tensor T and the energy flux j e . All of these assumptions/simplifications are common in kinetic-type theories of dilute polymeric fluids, and we also adopt them in the current contribution.
For further reference, we can therefore write the balance equations as
d ρ s d t + ρ s div x v = 0 ,
ρ s d v d t = div x T + ρ s b ,
ρ s d e d t = T : D div x j e ,
and we recall that the unknown fields ρ s and v are functions of the spatial position x and time t.
We note that, if one wanted to develop a model wherein the mass of the polymeric chains is not negligible, then one would need to use a variant of mixture theory, see, for example, Rajagopal and Tao [24] or Hutter and Jöhnk [25], and also Souček et al. [26] for a careful discussion.

3. Fokker–Planck Equation

We are in a position to describe the dynamics of the polymeric chains dispersed in the solvent. In what follows, we deal with a very simple setting, and we assume that the whole polymeric chain can be modelled as a single dumbbell (two beads connected with a spring), while more sophisticated models are typically necessary in order to obtain quantitative agreement between the model predictions and the experimental data (see the aforementioned monographs for a list of such models. Additional references regarding viscoelastic models can be also found in Vinogradov and Malkin [27], Larson [28] or Leonov and Prokunin [29]. A historical perspective is given in Tanner and Walters [30]). Note however that a generalisation to more complex polymeric chains/more sophisticated models is straightforward.
The dynamics of the individual polymeric chain are governed by a stochastic differential equation, which has the associated Fokker–Planck equation, see, for example, Le Bris and Lelievre [10] or Barrett and Süli [31] for details. If we consider dumbbells, their configuration is fully determined by the end-to-end vector q , and the Fokker–Planck equation governs the evolution of the configurational distribution function φ ( t , x , q ) in the configurational space D of all admissible end-to-end vectors. The physical meaning of the configurational distribution function is clear, the integral
D D φ ( t , x , q ) d q
gives the number of polymer chains at time t and spatial location x whose configuration vector lies in the domain D . The integral over the whole D gives the particle number density n p ( t , x ) ,
n p ( t , x ) = def D φ ( t , x , q ) d q ,
and the integral over a spatial domain
Ω D φ ( t , x , q ) d q d x ,
gives one the number of polymeric chains in the spatial domain Ω .
In the case of Hookean dumbbells, the configurational space D is the whole R 3 , since the stretch of the spring is not limited. In the case of the FENE dumbbell model (finitely extensible nonlinear elastic spring) the configurational space D is a ball in R 3 centred at the origin with a given radius, and the radius determines the maximum stretch of the spring.

3.1. Fokker–Planck Equation in the Case of Velocity Field with Nonzero Divergence

The Fokker–Planck equation has, in our case, the form
φ t + div x ( v φ + j φ , x ) + div q ( ( x v ) q φ + j φ , q ) = 0 ,
where j φ , x and j φ , q are unknown flux terms. Our task will be to identify these flux terms using the known energy storage mechanisms and entropy production mechanisms. We note that the spatially dependent term reads div x ( v φ ) and not v div x φ (this subtle difference does not matter in the case of incompressible fluids, but it is critical provided that we are working with a compressible fluid. See. for example, Degond and Liu [32] for a careful discussion of the Fokker–Planck equation). For further reference we also record that (7) can be rewritten as
d φ d t + φ div x v + div x j φ , x + div q ( x v ) q φ + j φ , q = 0 .

3.2. Boundary Condition in the Configurational Space

The boundary condition on D is the no-flux boundary condition
x v q φ + j φ , q n q D = 0 ,
where n q denotes the unit outward normal to the set of admissible end-to-end vectors D R 3 (in the case of Hookean dumbbells where D = R 3 the boundary condition is interpreted as a decay of the corresponding function at infinity). This choice of boundary condition has several implications regarding the evolution equation for polymer number density, see below, and regarding the behaviour of the configurational distribution function on the boundary of D, see Section 3.4.

3.3. Evolution Equation for Polymer Number Density

The evolution of the polymer number density n p introduced in (5) is governed by the partial differential equation
n p t + div x v n p + D j φ , x d q = 0 .
This equation is straightforward to obtain via integration of the Fokker–Plank equation with respect to the configurational variable q . The boundary term in the configurational space vanishes by virtue of (9). The last equation can be rewritten as
d n p d t + n p div x v = div x D j φ , x d q ,
where d d t denotes the standard material time derivative. Equation (11) shows that if j φ , x 0 and if the velocity field does not identically satisfy div x v = 0 , then, in general, we cannot expect that the polymer number density n p at the given material point is preserved. This makes the study of the dynamics of a compressible fluid with polymeric centre-of-mass diffusion substantially different from the simpler case of an incompressible fluid without polymeric centre-of-mass diffusion, see, for example, Barrett and Süli [33] and Feireisl et al. [34].

3.4. Force Potential

The force F in the spring connecting the beads is determined by a spherically symmetric potential
F = def q U 1 2 q q ref 2 ,
where q ref is a characteristic length scale for the spring. (Note the sign convention in (12) and the sign convention used later in (15).) For example, a popular choice for the force potential is the potential that is introduced by Warner [35] for the FENE dumbbell model. The potential is given by the formula
U FENE ( s ) = def b 2 log 1 2 s b .
We note that the choice (13) predicts infinite force as s b 2 . This means that an infinite force is necessary to maximally stretch the spring, which is the expected behaviour in view of the fact that we are dealing with a finitely extensible spring. The spring can not be extended beyond a certain length no matter what force is applied.
In what follows, we will not work with a specific force potential, we will only assume that the force potential U is composed of two parts U = U e + U η , θ , where
U e = U e 1 2 q q ref 2 ,
U η , θ = θ s θ ref U η 1 2 q q ref 2 ,
and where θ ref is a reference temperature, and q ref is a characteristic length scale. The potential U η , θ is the “entropic” part of the spring potential, and it is proportional to the ambient temperature θ s , while the potential U e is independent of the temperature (see, for example, Ericksen [36] for the rationale of the nomenclature, and the fact that it makes sense to consider the potential proportional to the temperature). If we are dealing with a finitely extensible spring, then the force tends to infinity as the spring reaches its maximum length (see for example the FENE potential introduced in (13)). This suggests that the value of the configurational distribution function φ on the boundary of D should be zero. We informally show that this property is in fact encoded in the no-flux boundary condition (9).
If we for a moment suppose that the Fokker–Planck Equation (7) reads
φ t + div x v φ k B θ s 2 ζ x φ + div q ( x v ) q φ 2 F ζ φ 2 k B θ s ζ q φ = 0 ,
that is if we consider specific formulae for the fluxes j φ , q and j φ , x , then (9) reduces to
( x v ) q φ 2 F ζ φ 2 k B θ s ζ q φ n q D = 0 .
Here, θ s denotes the solvent temperature, ζ denotes the hydrodynamic drag coefficient and k B denotes the Boltzmann constant. If D is a ball in R 3 , then q D = q max n q , and (16) reduces to
( x v ) q max n q φ 2 F ζ φ 2 k B θ s ζ q φ n q D = 0 .
Furthermore, in the case of spherically symmetric potential, the force F is parallel to q , and if we are dealing with a finitely extensible spring, then the potential is chosen in such a way that the magnitude of F approaches infinity as the vector q approaches the boundary of D. Having this piece of information in hand we inspect (17), and we see that if (17) holds, then the middle term 2 F ζ φ n q constrains the configurational distribution function φ to decay to zero faster than the growth of the force F as one approaches the boundary of D. Consequently, we can claim that (17) implies φ D = 0 .

3.5. Stationary Solution of the Fokker–Planck Equation in a Spatially Homogeneous State at a Given Temperature

If we assume for a moment that the fluxes in the general Fokker–Planck equation have been successfully identified, and have the form discussed in (15), then we can explicitly identify the stationary solution to the Fokker–Planck equation in the case of a spatially homogeneous temperature field and configurational distribution function field (partial derivatives with respect to time and the spatial variable x vanish). Under these assumptions, we see that the Fokker–Planck Equation (15) reduces to
div q 2 F ζ φ + 2 k B θ s ζ q φ = 0 .
If we assume that the force in the spring is given in terms of a potential, see (12), and if we substitute (12) into (18), then we get
div q q U 1 2 q q ref 2 φ + k B θ s q φ = 0 .
Using the formulae for the differential operators in a spherical coordinate system, we see that (19) reduces to a single ordinary differential equation
1 q 2 d d q q 2 d d q U 1 2 q q ref 2 φ + k B θ s d φ d q = 0 ,
where we have used the notation q = def q . We immediately see that the solution to (20) reads
φ = C e U 1 2 q q ref 2 k B θ s ,
where C is a constant. The constant is fixed by the condition n p ( t , x ) = D φ ( t , x , q ) d q , with φ now independent of t and x here, and hence the configurational distribution function that solves the stationary Fokker–Planck in the case of spatially homogeneous fields is given by the formula φ = M n p , θ s , where
M n p , θ s = def n p e U 1 2 q q ref 2 k B θ s D e U 1 2 q q ref 2 k B θ s d q ,
and n p is a positive constant. This configurational distribution function is the expected “equilibrium” configurational distribution function in the case when the dilute polymeric fluid occupies a thermodynamically isolated vessel. We return to this issue in Section 6.

4. Helmholtz Free Energy

Having discussed the Fokker–Planck equation that governs the evolution of the configurational distribution function, we have completed the system of equations for the unknown fields ρ s , v , θ s and φ . It remains to identify the fluxes, which are quantities that are unique for the given fluid (they provide a unique characterisation of the given material or given class of materials). The first step in this direction is the identification of energy storage mechanism, which is done by means of specific Helmholtz free energy.
We assume that the Helmholtz free energy ψ , ψ = J / kg , is given by the formula
ψ ( θ s , ρ s , φ ) = def ψ s ( θ s , ρ s ) + 1 ρ s D U e 1 2 q q ref 2 φ d q + k B θ s ρ s D U η 1 2 q q ref 2 k B θ ref φ + φ ln φ K d q ,
where K is a constant to be specified later (so far its only purpose is to provide a normalisation factor in the argument of the logarithm function, since we must have a dimensionless quantity under the logarithm. We also recall the physical units of the configurational distribution function φ , φ = 1 / m 6 , particle number density n p , n p = 1 / m 3 , potentials U e , U e = J , and U η , U η = J , and the Boltzmann constant k B , k B = J / K ). The inclusion of the potentials U e and U η into the formula for the Helmholtz free energy is a standard one known from elasticity theory, while the term φ ln φ is a well known contribution due to configurational entropy of the dumbbells. The solvent contribution ψ s ( θ s , ρ s ) to the specific Helmholtz free energy is left unspecified, since its particular form is not important for the ongoing discussion. We however give a particular formula for ψ s ( θ s , ρ s ) in Section 6, where it plays a fundamental role.
Having introduced the formula for the Helmholtz free energy, we can find other thermodynamic quantities of interest via differentiation. Namely, the formula for the specific entropy η reads
η ( θ s , ρ s , φ ) = ψ θ s = ψ s ( θ s , ρ s ) θ s k B ρ s D U η 1 2 q q ref 2 k B θ ref φ + φ ln φ K d q ,
and the specific internal energy e is given by the formula
e ( θ s , ρ s , φ ) = ψ + θ s η = ψ s ( θ s , ρ s ) θ s ψ s ( θ s , ρ s ) θ s + 1 ρ s D U e 1 2 q q ref 2 φ d q .

5. Constitutive Relations

Once we have specified the Helmholtz free energy we can proceed with the derivation of constitutive relations for the fluxes j e , j η , j φ , x , j φ , q , and the Cauchy stress tensor T . Our first objective is to identify the entropy flux j η and the entropy production ξ , such that the evolution equation for the specific entropy has the form
ρ s d η d t + div x j η = ξ .
The requirement on the nonnegativity of the entropy production then gives us a hint regarding the choice of constitutive relations for the fluxes j e , j η , j φ , x , j φ , q .
So far, we have identified an evolution equation for the internal energy (3c),
ρ s d e d t = T : D div x j e ,
and we know that the relation between the Helmholtz free energy and the internal energy and the entropy reads ψ ( θ s , ρ s , φ ) = e θ s η . Differentiation of this relation yields
ψ θ s d θ s d t + ψ ρ s d ρ s d t + ψ φ d φ d t = d e d t η d θ s d t θ s d η d t .
Making use of (27) and the fact that η = ψ θ s , we then see that (28) reduces to
ρ s θ s d η d t = T : D div x j e ρ s ψ ρ s d ρ s d t ρ s ψ φ d φ d t .
Now we need to manipulate (29) into the form (26). The lengthy but essential algebraic manipulations that allow us to achieve this objective are described in the following section, see Section 5.1. In particular, the desired final form of the entropy evolution equation is given in (51). The reader who is at the moment not interested in the necessary algebraic manipulations can go directly to Section 5.2, where we discuss the implications of the entropy evolution Equation (51) with respect to the choice of constitutive relations.

5.1. Evolution Equation for the Specific Entropy

Let us for a moment denote
F θ s , φ = def D U e 1 2 q q ref 2 φ + θ s θ ref U η 1 2 q q ref 2 φ + k B θ s φ ln φ K d q .
Using this notation the formula for the Helmholtz free energy reads
ψ ( θ s , ρ s , φ ) = ψ s ( θ s , ρ s ) + 1 ρ s F θ s , φ ,
and the evolution equation for the entropy (29) can be—with a slight abuse of notation—rewritten as
ρ s θ s d η d t = T : D div x j e ρ s ψ s ρ s d ρ s d t + F ρ s d ρ s d t F φ d φ d t .
Let us first consider the term F φ d φ d t . Making use of the Fokker–Planck Equation (8), we see that
F φ d φ d t = D U e 1 2 q q ref 2 + θ s θ ref U η 1 2 q q ref 2 + k B θ s ln φ K + 1 d φ d t d q = D U e 1 2 q q ref 2 + θ s θ ref U η 1 2 q q ref 2 + k B θ s ln φ K + 1 φ div x v + div x j φ , x + div q ( x v ) q φ + j φ , q d q = D U e 1 2 q q ref 2 + θ s θ ref U η 1 2 q q ref 2 + k B θ s ln φ K + 1 φ div x v + div x j φ , x d q A D U e 1 2 q q ref 2 + θ s θ ref U η 1 2 q q ref 2 + k B θ s ln φ K + 1 div q ( x v ) q φ + j φ , q d q B .
We will first focus on the term B that contains the divergence operator with respect to the q variable. We see that
D U e 1 2 q q ref 2 + θ s θ ref U η 1 2 q q ref 2 + k B θ s ln φ n p K + 1 div q ( x v ) q φ + j φ , q d q = D div q U e 1 2 q q ref 2 + θ s θ ref U η 1 2 q q ref 2 + k B θ s ln φ K + 1 ( x v ) q φ + j φ , q d q D q U e 1 2 q q ref 2 + θ s θ ref U η 1 2 q q ref 2 + k B θ s ln φ K + 1 ( x v ) q φ + j φ , q d q .
The first integral vanishes by virtue of the Stokes theorem and the boundary condition (9) on D . In the second term we first evaluate the gradient, and then we proceed, as follows
D q U e 1 2 q q ref 2 + θ s θ ref U η 1 2 q q ref 2 + k B θ s ln φ K + 1 ( x v ) q φ + j φ , q d q = D d U e d s s = 1 2 q q ref 2 + θ s θ ref d U η d s s = 1 2 q q ref 2 q q ref 2 + k B θ s q φ φ ( x v ) q φ + j φ , q d q = D d U e d s s = 1 2 q q ref 2 + θ s θ ref d U η d s s = 1 2 q q ref 2 q q ref q q ref φ : ( x v ) d q + D k B θ s q φ ( x v ) q d q + D 1 φ d U e d s s = 1 2 q q ref 2 + θ s θ ref d U η d s s = 1 2 q q ref 2 q q ref 2 φ + k B θ s q φ j φ , q d q = D d U e d s s = 1 2 q q ref 2 + θ s θ ref d U η d s s = 1 2 q q ref 2 q q ref q q ref φ : D d q + D k B θ s q φ ( x v ) q d q + D 1 φ d U e d s s = 1 2 q q ref 2 + θ s θ ref d U η d s s = 1 2 q q ref 2 q q ref 2 φ + k B θ s q φ j φ , q d q .
Furthermore, the second term on the right-hand side can be manipulated, as follows
D k B θ s q φ ( x v ) q d q = D div q k B θ s φ ( x v ) q d q D k B θ s φ x v : q q d q = D k B θ s φ div x v d q = k B θ s n p div x v ,
where we have again used the fact that the Stokes theorem and the boundary condition in configurational space (9) imply that the first term vanishes. As a matter of fact, we have used that φ D = 0 , which is a consequence of the particular choice of the flux j φ , q , see Section 3.4 for details. So far, we do not, however, have a formula for the flux j φ , q , and therefore we can not decide whether the term really vanishes. We will simply assume that the term vanishes, derive the formula for the flux j φ , q , and then we retrospectively check, that the derived formula for the flux j φ , q and boundary condition (9) are consistent with this assumption (this is indeed the case). We also note that the critical term can also be rewritten as
D div q k B θ s φ ( x v ) q d q = D k B θ s φ ( x v ) q n q d q = D k B θ s q max φ ( x v ) n q n q d q = D k B θ s q max φ n q n q d q : D ,
hence, if we wanted to keep it in the equations, we could easily do so. If we summarise the partial results, we therefore see that
B = D d U e d s s = 1 2 q q ref 2 + θ s θ ref d U η d s s = 1 2 q q ref 2 q q ref q q ref φ : D d q + k B θ s n p div x v D 1 φ d U e d s s = 1 2 q q ref 2 + θ s θ ref d U η d s s = 1 2 q q ref 2 q q ref 2 φ + k B θ s q φ j φ , q d q .
Next we focus on the term A in (33). We see that
D U e 1 2 q q ref 2 + θ s θ ref U η 1 2 q q ref 2 + k B θ s ln φ K + 1 φ div x v + div x j φ , x d q = D U e 1 2 q q ref 2 φ + θ s θ ref U η 1 2 q q ref 2 φ + k B θ s φ ln φ K + k B θ s φ div x v d q + D U e 1 2 q q ref 2 + θ s θ ref U η 1 2 q q ref 2 + k B θ s ln φ K + 1 div x j φ , x d q = F div x v + k B θ s n p div x v + D U e 1 2 q q ref 2 + θ s θ ref U η 1 2 q q ref 2 + k B θ s ln φ K + 1 div x j φ , x d q ,
and, hence
A = F div x v + k B θ s n p div x v + D U e 1 2 q q ref 2 + θ s θ ref U η 1 2 q q ref 2 + k B θ s ln φ K + 1 div x j φ , x d q .
Now, we can finally return to the formula (32) for the material time derivative of entropy, and we find that
ρ s θ s d η d t = T : D div x j e ρ s ψ s ρ s d ρ s d t + F ρ s d ρ s d t F φ d φ d t = T : D div x j e + ρ s 2 ψ s ρ s div x v F div x v + A + B ,
which reduces to
ρ s θ s d η d t = T : D div x j e + ρ s 2 ψ s ρ s div x v + 2 k B θ s n p div x v D d U e d s s = 1 2 q q ref 2 + θ s θ ref d U η d s s = 1 2 q q ref 2 q q ref q q ref φ : D d q D 1 φ d U e d s s = 1 2 q q ref 2 + θ s θ ref d U η d s s = 1 2 q q ref 2 q q ref 2 φ + k B θ s q φ j φ , q d q + D U e 1 2 q q ref 2 + θ s θ ref U η 1 2 q q ref 2 + k B θ s ln φ K + 1 div x j φ , x d q ,
where we have used the evolution equation for the solvent density (3a). We divide (42) by the solvent temperature, and manipulate the energy flux term div x j e in the usual manner. We also split the term containing div x j φ , x into a part with a temperature-independent coefficient and a part with a temperature-dependent coefficient, and after some manipulations we get
ρ s d η d t = 1 θ s T D d U e d s s = 1 2 q q ref 2 + θ s θ ref d U η d s s = 1 2 q q ref 2 q q ref q q ref φ : D + 1 θ s ρ s 2 ψ s ρ s + 2 k B θ s n p div x v 1 θ s D 1 φ d U e d s s = 1 2 q q ref 2 + θ s θ ref d U η d s s = 1 2 q q ref 2 q q ref 2 φ + k B θ s q φ j φ , q d q + D U e 1 2 q q ref 2 θ s div x j φ , x d q + D U η 1 2 q q ref 2 θ ref + k B ln φ K + 1 div x j φ , x d q div x j e θ s j e x θ s θ s 2 .
We introduce the notation
P = def D d U e d s s = 1 2 q q ref 2 + θ s θ ref d U η d s s = 1 2 q q ref 2 q q ref q q ref φ d q ,
and we split the Cauchy stress tensor T into the mean normal stress m and the traceless part T δ ,
T = m I + T δ ,
where m = def 1 3 Tr T and T δ = def T 1 3 Tr T I (this standard manipulation allows us to treat the volume-changing and volume-preserving deformations separately, see for example Málek and Průša [21] for details). We also do the same for the tensor P , and we rewrite (43) as
ρ s d η d t = 1 θ s T δ P δ : D δ + m 1 3 Tr P + ρ s 2 ψ s ρ s + 2 k B θ s n p div x v 1 θ s D 1 φ d U e d s s = 1 2 q q ref 2 + θ s θ ref d U η d s s = 1 2 q q ref 2 q q ref 2 φ + k B θ s q φ j φ , q d q + D U e 1 2 q q ref 2 θ s div x j φ , x d q + D U η 1 2 q q ref 2 θ ref + k B ln φ K + 1 div x j φ , x d q div x j e θ s j e x θ s θ s 2 .
Now, we are ready to manipulate the part with the temperature independent coefficient
D U η 1 2 q q ref 2 θ ref + k B ln φ K + 1 div x j φ , x d q = div x D U η 1 2 q q ref 2 θ ref + k B ln φ K + 1 j φ , x d q D x U η 1 2 q q ref 2 θ ref + k B ln φ K + 1 j φ , x d q .
Substituting (47) back into (43) yields
ρ s d η d t = 1 θ s T δ P δ : D δ + m 1 3 Tr P + ρ s 2 ψ s ρ s + 2 k B θ s n p div x v 1 θ s D 1 φ d U e d s s = 1 2 q q ref 2 + θ s θ ref d U η d s s = 1 2 q q ref 2 q q ref 2 φ + k B θ s q φ j φ , q d q + D U e 1 2 q q ref 2 θ s div x j φ , x d q D k B φ x φ j φ , x d q div x j e θ s D U η 1 2 q q ref 2 θ ref + k B ln φ K + 1 j φ , x d q j e x θ s θ s 2 .
We rewrite the last equation as
ρ s d η d t + div x j e θ s D U η 1 2 q q ref 2 θ ref + k B ln φ K + 1 j φ , x d q = 1 θ s T δ P δ : D δ + m 1 3 Tr P + ρ s 2 ψ s ρ s + 2 k B θ s n p div x v 1 θ s D 1 φ d U e d s s = 1 2 q q ref 2 + θ s θ ref d U η d s s = 1 2 q q ref 2 q q ref 2 φ + k B θ s q φ j φ , q d q + D U e 1 2 q q ref 2 θ s div x j φ , x d q D k B φ x φ j φ , x d q j e x θ s θ s 2 ,
which helps us to identify the fluxes and entropy production terms. It remains to manipulate the term with the potential U e . We see that
D U e 1 2 q q ref 2 θ s div x j φ , x d q = div x D U e 1 2 q q ref 2 θ s j φ , x d q + x θ s D U e 1 2 q q ref 2 j φ , x d q θ s 2 .
This manipulation implies that we can rewrite (49) in the final form as
ρ s d η d t + div x j e θ s D U η 1 2 q q ref 2 θ ref + k B ln φ K + 1 + U e 1 2 q q ref 2 θ s j φ , x d q = 1 θ s T δ P δ : D δ + m 1 3 Tr P + ρ s 2 ψ s ρ s + 2 k B θ s n p div x v 1 θ s D 1 φ d U e d s s = 1 2 q q ref 2 + θ s θ ref d U η d s s = 1 2 q q ref 2 q q ref 2 φ + k B θ s q φ j φ , q d q D k B φ x φ j φ , x d q j e D U e 1 2 q q ref 2 j φ , x d q x θ s θ s 2 .

5.2. Entropy Production and Constitutive Relations

Now, we are in a position to choose j η , j e , j φ , x , j φ , q , m and T δ in such a way that the right-hand side is nonnegative. We start with the flux j φ , x in physical space, the nonnegativity of the term D k B φ x φ j φ , x d q is granted if we make, for example, the simple choice
j φ , x = def k B θ s 2 ζ x φ ,
where ζ denotes the hydrodynamic drag coefficient, ζ = N m s (the presence of the hydrodynamic drag coefficient is motivated by the insight into the microscopic dynamics of the polymeric chains, see for example Barrett and Süli [31] for details). The choice (52) in turn implies that we have to fix the energy flux as
j e = def κ x θ s + D U e 1 2 q q ref 2 j φ , x d q = κ x θ s k B θ s 2 ζ D U e 1 2 q q ref 2 x φ d q .
We note that (53) is the same energy flux as in Öttinger and Grmela [18]. We also note that only the “energetic” part U e of the potential U enters the formula for the energy flux j e , which is an expected result.
Next, we fix the constitutive relation for the flux j φ , q in the configurational space; we set
j φ , q = def 2 ζ d U e d s s = 1 2 q q ref 2 + θ s θ ref d U η d s s = 1 2 q q ref 2 q q ref 2 φ + k B θ s q φ ,
where the coefficient 2 ζ is again motivated by the microscopic dynamics of the polymeric chains. In particular, we want to keep the link between the Fokker–Planck equation and the associated stochastic process, see, for example, Barrett and Süli [31] for details. If we recall the relation between the force and the corresponding potential (12), we see that (54), in fact, reads
j φ , q = 2 ζ F φ + k B θ s q φ ,
where
F = def q U e 1 2 q q ref 2 + θ s θ ref U η 1 2 q q ref 2 ,
which, together with (52), gives upon substitution into (7) the Fokker–Planck equation in the form
φ t + div x v φ k B θ s 2 ζ x φ + div q ( x v ) q φ 2 F ζ φ 2 k B θ s ζ q φ = 0 .
This is the specific Fokker–Planck equation that we have used in Section 3.
From the left-hand side of (51) we can infer the entropy flux j η . If we use the formula for j e , we see that the term under the div x operator reads
j η = def κ x θ s θ s D U η 1 2 q q ref 2 θ ref + k B ln φ K + 1 j φ , x d q ,
where j φ , x is given by (52). We note that the term with the “energetic” part U e of the potential U has cancelled, and only the entropic part U η enters the entropic flux, which is an expected result.
Finally, we fix the constitutive relation for the Cauchy stress tensor, namely for the mean normal stress m and the deviatoric part T δ . We need a choice that makes the term T δ P δ : D δ + m 1 3 Tr P + ρ s 2 ψ s ρ s + 2 k B θ s n p div x v on the right-hand side of (51) nonnegative. If we want to model the solvent as a classical compressible Navier–Stokes fluid, then we set
T δ = def 2 ν D δ + P δ ,
m = def 1 3 Tr P ρ s 2 ψ s ρ s 2 k B θ s n p + λ ˜ div x v ,
where ν and λ ˜ are the positive constants referred to as the shear viscosity and bulk viscosity. If we use (59), then the full formula for the Cauchy stress tensor reads
T = ρ s 2 ψ s ρ s I + λ div x v I + 2 ν D 2 k B θ s n p I + P ,
where λ ˜ = 3 λ + 2 ν 3 . The first three terms on the right-hand side of (60) are the familiar ones. In particular, we recognise the formula for the thermodynamic pressure
p th , s ( θ s , ρ s ) = def ρ s 2 ψ s ρ s ;
hence, we can rewrite (60) as
T = p th , s I + λ div x v I + 2 ν D 2 k B θ s n p I + P .
The polymeric contribution T polymer to the Cauchy stress tensor T is given by the last two terms in (62). If we recall the definition of the tensor P , see (44), we see that
T polymer = 2 k B θ s n p I + D d U e d s s = 1 2 q q ref 2 + θ s θ ref d U η d s s = 1 2 q q ref 2 q q ref q q ref φ d q .
This is the Kramers expression for the polymeric contribution to the Cauchy stress tensor. We note that the Kramers expression has not been assumed a priori: it is a consequence of the choice of energy storage and entropy production mechanisms. The same formula also follows from the GENERIC formalism, see Öttinger and Grmela [18], where it is also a consequence of the proposed modelling approach and it is not assumed a priori.

5.3. Temperature Evolution Equation

A cautious reader might notice that the energy flux j e and the entropy flux j η depend on the potentials U e and U η and on the normalisation constant K. However, if one shifts the value of a potential by a constant, then this operation should have no impact on the dynamics of the system. The same should hold for the constant K. Its particular choice should have no impact on the dynamics.
We show that this is indeed true. First, we observe that the shift of the value of the potentials and the choice of K have no impact on the balance of mass (3a) (evolution equation for the solvent density ρ s ) and on the balance of linear momentum (78b) (evolution equation for the velocity field v ). Indeed, the Cauchy stress tensor, see (60), depends on the gradients of the potentials. Next we observe that the Fokker–Planck equation (evolution equation for φ ) in the form (57) also contains only gradients of the potentials. It remains to check the evolution equation for the last quantity of interest, namely for the temperature.
In order to do so, we first need to explicitly formulate the evolution equation for the temperature. We start with the evolution equation for the entropy in the form (29), which is
ρ s θ s d η d t = T : D div x j e ρ s ψ ρ s d ρ s d t ρ s ψ φ d φ d t .
Because the entropy is given as the derivative of the Helmholtz free energy with respect to the temperature, η = ψ θ s , we can use the chain rule on the left-hand side of (64), and with a slight abuse of notation we get
ρ s θ s d η d t = ρ s θ s 2 ψ θ s 2 d θ s d t 2 ψ θ s ρ s d ρ s d t 2 ψ θ s φ d φ d t .
Making use of (65) in (64) yields
ρ s θ s 2 ψ θ s 2 d θ s d t = ρ s θ s 2 ψ θ s ρ s d ρ s d t + ρ s θ s 2 ψ θ s φ d φ d t + T : D div x j e ρ s ψ ρ s d ρ s d t ρ s ψ φ d φ d t .
On the left-hand side, we have obtained the material time derivative of the temperature field θ s , and we can explicitly evaluate all terms on the right-hand side. We take into account the particular structure of the Helmholtz free energy (23), which is
ψ ( θ s , ρ s , φ ) = def ψ s ( θ s , ρ s ) + 1 ρ s D U e 1 2 q q ref 2 φ d q + k B θ s ρ s D U η 1 2 q q ref 2 k B θ ref φ + φ ln φ K d q ψ p
and we arrive at
ρ s θ s 2 ψ s θ s 2 d θ s d t = ρ s θ s 2 ψ s θ s ρ s ψ s ρ s d ρ s d t + T : D div x j e + ρ s θ s 2 ψ p θ s ρ s d ρ s d t + ρ s θ s 2 ψ p θ s φ d φ d t ρ s ψ p ρ s d ρ s d t ρ s ψ p φ d φ d t .
By virtue of the linearity of the polymeric part ψ p with respect to temperature, we get
ρ s θ s 2 ψ s θ s 2 d θ s d t = ρ s θ s 2 ψ s θ s ρ s ψ s ρ s d ρ s d t + T : D div x j e + 1 ρ s D U e 1 2 q q ref 2 φ d q d ρ s d t D U e 1 2 q q ref 2 d φ d t d q .
Next, we recall that the expression on the left-hand side of (69) is the classical definition of the specific heat at constant volume c V , s for the solvent only,
c V , s ( θ s , ρ s ) = def θ s 2 ψ s θ s 2 ,
we use the definition of the thermodynamic pressure p th , s for the solvent, see (61), the balance of mass (3a), and we get
ρ s c V , s d θ s d t = θ s p th , s θ s div x v + p th , s div x v + T : D div x j e D U e 1 2 q q ref 2 φ d q div x v D U e 1 2 q q ref 2 d φ d t d q .
Now, we are in a position to use the formula for the energy flux j e , see (53), and the Fokker–Planck Equation (8) for the configurational distribution function φ , which yields
ρ s c V , s d θ s d t = θ s p th , s θ s div x v + p th , s div x v + T : D + div x κ x θ s + D U e 1 2 q q ref 2 div q ( x v ) q φ + j φ , q d q .
Let us now focus on the last term. Integration by parts in the last term together with the boundary condition (9) and the constitutive relation (54) for j φ , q then gives us
D U e 1 2 q q ref 2 div q ( x v ) q φ + j φ , q d q = D div q U e 1 2 q q ref 2 ( x v ) q φ + j φ , q d q D d U e d s s = 1 2 q q ref 2 q q ref 2 ( x v ) q φ + j φ , q d q = D d U e d s s = 1 2 q q ref 2 q q ref 2 ( x v ) q φ + j φ , q d q = D d U e d s s = 1 2 q q ref 2 q q ref q q ref φ d q : D D d U e d s s = 1 2 q q ref 2 q q ref 2 j φ , q d q = D d U e d s s = 1 2 q q ref 2 q q ref q q ref φ d q : D + 2 ζ D d U e d s s = 1 2 q q ref 2 q q ref 2 d U e d s s = 1 2 q q ref 2 + θ s θ ref d U η d s s = 1 2 q q ref 2 q q ref 2 φ + k B θ s q φ d q = D d U e d s s = 1 2 q q ref 2 q q ref q q ref φ d q : D + 2 ζ D d U e d s s = 1 2 q q ref 2 q q ref 2 d U e d s s = 1 2 q q ref 2 + θ s θ ref d U η d s s = 1 2 q q ref 2 q q ref 2 φ d q + 2 k B θ s ζ D d U e d s s = 1 2 q q ref 2 q q ref 2 q φ d q .
If we want to, we can also integrate by parts in the last term on the right-hand side, and we get
2 k B θ s ζ D d U e d s s = 1 2 q q ref 2 q q ref 2 q φ d q = 2 k B θ s ζ D div q d U e d s s = 1 2 q q ref 2 q q ref 2 φ d q 2 k B θ s ζ D div q d U e d s s = 1 2 q q ref 2 q q ref 2 φ d q = 2 k B θ s ζ D div q q U e 1 2 q q ref 2 φ d q = 2 k B θ s ζ D Δ q U e 1 2 q q ref 2 φ d q ,
where we have again exploited the fact that φ D = 0 . Now, we substitute back into (72), and we arrive at
ρ s c V , s d θ s d t = θ s p th , s θ s div x v + p th , s div x v + T : D + div x κ x θ s D d U e d s s = 1 2 q q ref 2 q q ref q q ref φ d q : D + 2 ζ D d U e d s s = 1 2 q q ref 2 q q ref 2 d U e d s s = 1 2 q q ref 2 + θ s θ ref d U η d s s = 1 2 q q ref 2 q q ref 2 φ d q 2 k B θ s ζ D Δ q U e 1 2 q q ref 2 φ d q .
Finally, we use the formula for the Cauchy stress tensor (62), which is
T = p th , s I + λ div x v I + 2 ν D 2 k B θ s n p I + D d U e d s s = 1 2 q q ref 2 + θ s θ ref d U η d s s = 1 2 q q ref 2 q q ref q q ref φ d q ,
which, upon substitution into (75), yields
ρ s c V , s d θ s d t = θ s p th , s θ s div x v + div x κ x θ s + λ div x v 2 + 2 ν D : D 2 k B θ s n p div x v + D θ s θ ref d U η d s s = 1 2 q q ref 2 q q ref q q ref φ d q : D + 2 ζ D d U e d s s = 1 2 q q ref 2 q q ref 2 d U e d s s = 1 2 q q ref 2 + θ s θ ref d U η d s s = 1 2 q q ref 2 q q ref 2 φ d q 2 k B θ s ζ D Δ q U e 1 2 q q ref 2 φ d q .
This is the sought evolution equation for the temperature field. A brief inspection of (77) reveals that the evolution equation for the temperature is not affected by shifts of the values of the potentials U e and U η and the value of scaling constant K.

5.4. Summary

The system of governing equations for the solvent density ρ s ( t , x ) , velocity v ( t , x ) , solvent temperature θ s ( t , x ) , and configurational distribution function φ ( t , x , q ) read, as follows:
d ρ s d t + ρ s div x v = 0 ,
ρ s d v d t = div x T + ρ s b ,
φ t + div x v φ k B θ s 2 ζ x φ + div q ( x v ) q φ 2 F ζ φ 2 k B θ s ζ q φ = 0 ,
where the spring force F is given via the potentials U e and U η by the formula
F = q U e 1 2 q q ref 2 + θ s θ ref U η 1 2 q q ref 2 = d U e d s s = 1 2 q q ref 2 + θ s θ ref d U η d s s = 1 2 q q ref 2 q q ref 2 ,
and the Cauchy stress tensor T is given by the formula
T = p th , s I + λ div x v I + 2 ν D 2 k B θ s n p I + D F q φ d q .
The temperature evolution equation reads
ρ s c V , s d θ s d t = θ s p th , s θ s div x v + div x κ x θ s + λ div x v 2 + 2 ν D : D 2 k B θ s n p div x v + D q θ s θ ref U η q φ d q : D + 2 ζ D q U e q U e + θ s θ ref U η d q 2 k B θ s ζ D Δ q U e φ d q ,
and the thermodynamic pressure p th , s and the specific heat at constant volume c V , s are given in terms of the derivative of the solvent specific Helmholtz free energy ψ s , as follows:
p th , s ( θ s , ρ s ) = ρ s 2 ψ s ρ s ,
c V , s ( θ s , ρ s ) = θ s 2 ψ s θ s 2 ,
We recall that the solvent density ρ s is in the given model identified with the density of the whole solvent-polymeric chains mixture. We note that if the spring potential only contains the entropic part, that is if U e = 0 , then the governing equations simplify considerably. This system of governing equations is essentially the same as that reported by Öttinger and Grmela [18], albeit the model has been derived in a different manner. In fact our approach is similar to that used in the derivation of macroscopic viscoelastic rate-type models, see especially Rajagopal and Srinivasa [37], Málek et al. [38], Hron et al. [17], Málek et al. [39], Málek et al. [40], and Dostalík et al. [41]. For a coupling between microscopic and macroscopic models, see also Souček et al. [42].
The boundary condition in the configurational space D reads
( x v ) q φ 2 F ζ φ 2 k B θ s ζ q φ n q D = 0 ,
see discussion in Section 3.2. Because we are working with a model with centre-of-mass diffusion, we also need a boundary condition for the configurational distribution function in the physical space. If we consider a vessel Ω with an impermeable rigid wall, we set
j φ , x n x Ω = 0 ,
which, by virtue of (52), means that x φ n x Ω = 0 . Regarding the remaining boundary conditions in physical space, we can opt for any standard boundary conditions. Because we have access to the velocity field v as well as to the Cauchy stress tensor T , we can choose from plethora of boundary conditions used for flows of polymeric fluids, see, for example, Hatzikiriakos [43] or Málek and Průša [21].

6. Stability

Following Coleman and Greenberg [44] and Coleman [19], see also Gurtin [45,46], Grmela and Öttinger [47] and Bulíček et al. [20] for further discussion, we can exploit the thermodynamic basis of the derived model in a rudimentary stability analysis of thermo-mechanical processes described by the corresponding governing equations. In particular, it is straightforward to investigate the finite amplitude (nonlinear) stability of the spatially homogeneous stationary states in a thermodynamically isolated container.
The spatially homogeneous equilibrium stationary state θ s , eq , ρ s , eq , φ eq and v eq is, in our case, given as
θ s , eq = θ s ^ ,
ρ s , eq = ρ s ^ ,
φ eq = M n p ^ , θ s ^ ,
v eq = 0 .
The symbol n p ^ is the spatially homogeneous polymer number density, which is n p ^ = N Ω , where N is the total number of polymeric chains in the container and Ω is the volume of the container, the solvent density ρ s ^ is the spatially homogeneous density, which is ρ s ^ = M Ω , where M is the total mass of the solvent in the container. The symbol θ s ^ denotes a spatially homogeneous temperature field, and M n p ^ , θ s ^ stands for the equilibrium configurational distribution function (22) discussed in Section 3.5. Clearly, the quadruple (81) is a solution to the governing equations subject to boundary conditions
v Ω = 0 ,
x φ n x Ω = 0 ,
x θ s n x Ω = 0 ,
that guarantee that the energy flux j e , j η , and j φ , x , see (53), (58), and (52), vanish on the container wall. We also have the standard boundary condition in the configurational space (9), which is
x v q φ + j φ , q n q D = 0 ,
where j φ , q is given by (54). If we denote
S θ ^ = def Ω ρ s θ s ^ η d v ,
S = def Ω ρ s η d v ,
E tot = def Ω ρ s e + 1 2 ρ s v 2 d v ,
then using the governing equations it is straightforward to see that, since all of the fluxes through the vessel wall vanish, the net total energy E tot is constant and the net total entropy S is a nondecreasing function in time. Moreover, the total mass of the solvent and the number of polymeric chains are also preserved in time.
Consequently, the following functional is a reasonable candidate for a Lyapunov type functional for the analysis of the stability of the stationary spatially homogeneous state with temperature θ s ^ ,
V meq = S θ s ^ E tot E tot ^ ,
see, for example, Bulíček et al. [20] for details, where E tot ^ denotes the net total energy at the spatially homogeneous stationary state (we use the nomenclature Lyapunov type functional, since we will only show the decay along trajectories and non-negativity of the functional everywhere except at the equilibrium. We will not investigate the link between the proposed functional and a suitable norm in the given state space, as this is beyond the scope of the current contribution). The rationale for the choice of such a functional is apparent from the following manipulation
d V meq d t = d d t S θ ^ E tot E tot ^ = θ ^ d S d t = θ ^ Ω ξ d v ,
where ξ denotes the entropy production, which is a nonnegative quantity that vanishes at the spatially homogeneous stationary state (81). In the current case, we are forced however to modify the functional by explicitly adding zero terms
V meq = S θ s ^ E tot E tot ^ + λ 2 Ω ρ s ρ s ^ d v + λ 3 Ω n p n p ^ d v ,
which help us to articulate the constraints of mass conservation and the conservation of total number of polymeric chains. While the values of Lagrange multipliers can be found a priori, see, for example, Coleman [19] for a discussion in the case of a compressible fluid, we shall not proceed in that direction (in fact we have already tacitly identified one of the multipliers: the equilibrium temperature θ s ^ is the multiplier that enforces the conservation of net total energy). We shall find instead the multipliers on-the-fly. Indeed, from the pragmatic point of view, we will need the additional zero terms (after the integration) to show the nonnegativity of the integrands in V meq .
Finally, we note that, in order to ensure that the functional vanishes at the equilibrium, we need to shift the functional by a constant value. We will denote the shifted functional by the same symbol, since this shift has no impact on the time evolution of the functional: it is just a matter of normalisation.

6.1. Outline of the Construction of the Lyapunov Like Functional

We see that the explicit formula for V meq reads
V meq = Ω g θ s , ρ s , φ θ s ^ , ρ s ^ , φ ^ + 1 2 ρ s v 2 d v
where
g θ s , ρ s , φ θ s ^ , ρ s ^ , φ ^ = def ρ s θ s ^ η ( θ s , ρ s , φ ) ρ s e ( θ s , ρ s , φ ) e ( θ s ^ , ρ s ^ , φ ^ ) .
By virtue of (24) and (25), we can further split the internal energy and the entropy into their respective solvent part and polymeric part:
e θ s , ρ s , φ = e s θ s , ρ s + e p θ s , ρ s , φ ,
η θ s , ρ s , φ = η s θ s , ρ s + η p θ s , ρ s , φ ,
where
e p θ s , ρ s , φ = def 1 ρ s D U e 1 2 q q ref 2 φ d q ,
η p θ s , ρ s , φ = def k B ρ s D U η 1 2 q q ref 2 k B θ ref φ + φ ln φ K d q .
Next, we split the function g into a solvent and a polymeric part as well
g θ s , ρ s , φ θ s ^ , ρ s ^ , φ ^ = g s θ s , ρ s θ s ^ , ρ s ^ , + g p θ s , ρ s , φ θ s ^ , ρ s ^ , φ ^ ,
where
g s θ s , ρ s θ s ^ , ρ s ^ = def ρ s θ s ^ η s ( θ s , ρ s ) ρ s e s ( θ s , ρ s ) e s ( θ s ^ , ρ s ^ ) ,
g p θ s , ρ s , φ θ s ^ , ρ s ^ , φ ^ = def ρ s θ s ^ η p ( θ s , ρ s , φ ) ρ s e p ( θ s , ρ s , φ ) e p ( θ s ^ , ρ s ^ , φ ^ ) ,
and, finally, we also split the whole functional, which is we write
V meq = Ω 1 2 ρ s v 2 d v + Ω g s θ s , ρ s θ s ^ , ρ s ^ d v V meq , s + Ω g p θ s , ρ s , φ θ s ^ , ρ s ^ , φ ^ d v V meq , p .
We shall first find an explicit formula for the polymeric contribution, see Section 6.2, and then we introduce an explicit equation of state for the solvent, and we find the formula for the solvent contribution, see Section 6.3.

6.2. Polymeric Part

First, we choose the characteristic temperature θ ref in the polymeric part of the Helmholtz free energy (23) as θ s ^ , which is we set
θ ref = def θ s ^ .
We are using the fact that the temperature field θ s ^ at equilibrium is constant in space and time. Straightforward substitution of (90) into (92b) then yields
g p θ s , ρ s , φ θ s ^ , ρ s ^ , φ ^ = k B θ s ^ D U η 1 2 q q ref 2 k B θ s ^ φ + φ ln φ K d q + D U e 1 2 q q ref 2 φ d q D U e 1 2 q q ref 2 φ ^ d q .
Now, we can combine the terms with the potentials to deduce that
g p θ s , ρ s , φ θ s ^ , ρ s ^ , φ ^ = k B θ s ^ D φ ln φ K e U η 1 2 q q ref 2 + U e 1 2 q q ref 2 k B θ s ^ d q D U e 1 2 q q ref 2 φ ^ d q .
Because the choice of the normalisation factor K has no impact on the dynamics of the system, we choose it to serve our objective. We fix the normalisation factor K as
K = def n p ^ D e U η 1 2 q q ref 2 + U e 1 2 q q ref 2 k B θ s ^ d q ,
Using the notation introduced in (22), we see that (96) with the constant K chosen as in (97) reduces to
g p θ s , ρ s , φ θ s ^ , ρ s ^ , φ ^ = k B θ s ^ D φ ln φ M n p ^ , θ s ^ d q D U e 1 2 q q ref 2 φ ^ d q .
The last term is a constant in time and space, hence we can safely ignore it. It only shifts the value of g p and consequently of the functional V meq . If we ignore it, we get in the end a functional that vanishes at equilibrium, which is a desired property. We see that
1 M n p ^ , θ s ^ = 1 M n p , θ s ^ n p n p ^ ,
and by virtue of the normalisation constraint (5), it follows that
k B θ s ^ D φ ln φ M n p ^ , θ s ^ d q = k B θ s ^ D φ ln φ M n p , θ s ^ n p n p ^ d q = k B θ s ^ D φ ln φ M n p , θ s ^ d q + k B θ s ^ D φ ln n p n p ^ d q = k B θ s ^ D M n p , θ s ^ φ M n p , θ s ^ ln φ M n p , θ s ^ φ M n p , θ s ^ + 1 d q + k B θ s ^ n p ln n p n p ^ .
It is straightforward to check that the function in the square bracket, f ( x ) = def x ln x x + 1 , is for x ( 0 , + ) nonnegative and vanishes if and only if x = 1 , which is if
φ M n p , θ s ^ = 1 .
Therefore, the corresponding term is related to the equilibration of the configurational distribution function φ to M n p , θ s ^ . We note that M n p , θ s ^ depends on the local particle number density n p and not on the equilibrium particle number density n p ^ . Consequently, from this term only, we cannot monitor the progress of n p to n p ^ . This piece of information must be inferred from elsewhere.
Now, we are in a position to construct the functional V meq , p . In order to guarantee the nonnegativity of the functional we need to add a term that reflects the additional constraints imposed on the closed system, namely, the constraint on the constant total number of polymeric chains. See the discussion in the introductory part of this section. We set
V meq , p = Ω g p θ s , ρ s , φ θ s ^ , ρ s ^ , φ ^ d v k B θ s ^ Ω n p n p ^ d v ,
which, by virtue of (98) and (100), yields
V meq , p = k B θ s ^ Ω D M n p , θ s ^ φ M n p , θ s ^ ln φ M n p , θ s ^ φ M n p , θ s ^ + 1 d q d v + k B θ s ^ Ω n p ^ n p n p ^ ln n p n p ^ n p n p ^ + 1 d v .
A brief inspection of the last term reveals that the term is nonnegative, and that it vanishes if and only if
n p n p ^ = 1
everywhere in the domain Ω . Consequently, we are done with the polymeric contribution to the functional V meq , p and we can proceed to the solvent contribution.

6.3. Solvent Part: Noble–Abel Stiffened-Gas Equation of State

Regarding the solvent part, we need to choose an equation of state. In the current contribution, we opt for the Noble–Abel stiffened-gas equation of state (NASG), see Le Métayer and Saurel [48], which is a recently formulated equation of state that is easy to deal with analytically, and that can describe, by a suitable choice of parameters, both gaseous and liquid substances (see also Chiapolino and Saurel [49] for further elaboration of NASG). Let us follow Le Métayer and Saurel [48], and let us first recall various formulae of the Noble–Abel stiffened-gas equations of state. The quantities of interest in engineering practice are given by the formulae
p th , s ( θ s , ρ s ) = c V , s γ 1 ρ s θ s 1 b ρ s p ,
e s ( p th , s , ρ s ) = p th , s + γ p p th , s + p c V , s θ s + q ,
η s ( p th , s , ρ s ) = c V , s ln θ s θ s , ref γ p th , s + p p ref 1 γ + q ,
where p th , s , θ s , ρ s , e s , and η s are, respectively, the thermodynamic pressure, the temperature, the density, the specific internal energy, and the specific entropy. The symbols c V , s , γ , b, p , q, and q denote constant material parameters that are specific for the given fluid. The symbols θ s , ref and p ref denote constant reference temperature and constant reference pressure. These quantities can be chosen at will, as they serve only for normalisation purposes. Combining the Formulae (105), we obtain
e s θ s , ρ s = c V , s θ s + 1 ρ s b p + q ,
η s θ s , ρ s = c V , s ln θ s θ s , ref + c V , s 1 γ ln c V , s γ 1 ρ s 1 b ρ s θ s , ref p ref + q .
Finally, using (106), we obtain the specific Helmholtz free energy ψ s defined by the thermodynamic relation ψ s ( θ s , ρ s ) = e s ( θ s , ρ s ) θ s η s ( θ s , ρ s ) . The explicit formula for the specific Helmholtz free energy as a function of its natural variables reads
ψ s ( θ s , ρ s ) = c V , s θ s ln θ s θ s , ref 1 c V , s θ s 1 γ ln c V , s γ 1 ρ s 1 b ρ s θ s , ref p ref q θ s + 1 ρ s b p + q ,
and using this formula in (23) gives us a complete thermodynamic characterisation of the given fluid.
We are now in a position to construct the functional V meq , s , see (93). As in the previous section, we first choose the characteristic temperature of the solvent θ s , ref as θ s ^ , which is we set
θ s , ref = def θ s ^ .
We are using the fact that the temperature field θ s ^ at equilibrium is constant in space and time. With this choice of the reference temperature field θ s , ref , we use the formulae for the specific internal energy e s , and the specific entropy η s , see (106), and we get
V meq , s = Ω ρ s c V , s θ s ^ ln θ s θ s ^ + ρ s c V , s 1 γ θ s ^ ln c V , s γ 1 ρ s 1 b ρ s θ s ^ p ref + ρ s q θ s ^ d v + Ω ρ s c V , s θ s + 1 b ρ s p + ρ s q d v Ω ρ s ^ c V , s θ s ^ + 1 b ρ s ^ p + ρ s ^ q d v .
Now, we need to manipulate the right-hand side of (109), and we need to show that it is nonnegative. This can be done, as follows. First, we exploit conservation of mass, which implies that Ω ρ s ρ s ^ d v = 0 . See also the discussion in the introductory part, in particular the discussion following (86). This observation allows for us to write Ω ρ s c V , s θ s ^ d v = Ω ρ s ^ c V , s θ s ^ d v . Second, we fix the reference pressure as
p ref = def ρ s ^ γ 1 c V , s θ s ^ 1 b ρ s ^ e q c V , s γ 1 ,
and regroup the terms of (109) arriving at
V meq , s = Ω ρ s c V , s θ s ^ θ s θ s ^ 1 ln θ s θ s ^ d v + Ω ρ s c V , s γ 1 θ s ^ ln ρ s ρ s ^ 1 b ρ s ^ 1 b ρ s d v + Ω ρ s ρ s ^ q b p d v .
We note that the choice (110) leads to vanishing entropy at the equilibrium state. In this sense, the choice (110) is tantamount to the convenient fix of the entropy value at the equilibrium. Next we again exploit conservation of mass, which allows us to eliminate the last term on the right-hand side of (111), and yields
V meq , s = Ω ρ s c V , s θ s ^ θ s θ s ^ 1 ln θ s θ s ^ d v + Ω ρ s c V , s γ 1 θ s ^ ln ρ s ρ s ^ 1 b ρ s ^ 1 b ρ s d v .
Moreover, by adding the zero term
Ω c V , s γ 1 θ s ^ 1 b ρ s ^ ρ s ρ s ^ d v ,
to the right-hand side of (112), we arrive at
V meq , s = Ω ρ s c V , s θ s ^ θ s θ s ^ 1 ln θ s θ s ^ d v + Ω c V , s ( γ 1 ) θ s ^ ρ s ln ρ s ρ s ^ 1 b ρ s ^ 1 b ρ s ρ s ρ s ^ 1 b ρ s ^ d v .
Because the function f ( x ) = def x 1 ln x , is nonnegative for x 0 , and vanishes if and only if x = 1 , we obtain the desired nonnegativity of the first term on the right-hand side of (114).
Finally, we denote
g ρ s ^ , b ( ρ s ) = def ρ s ln ρ s ρ s ^ 1 b ρ s ^ 1 b ρ s ρ s ρ s ^ 1 b ρ s ^ ,
where ρ s 0 , 1 b . Clearly, g ρ s ^ , b ( ρ s ^ ) = 0 . Moreover, by calculating the first and the second derivatives of the function g ρ s ^ , b one obtains that, irrespective of the values of b > 0 and ρ s ^ 0 , 1 b , the function g ρ s ^ , b is nonnegative and vanishes if and only if ρ s = ρ s ^ . The second term on the right-hand side of (114) is thus nonnegative as well. We also note that the function g ρ s ^ , b ( ρ s ) blows-up as ρ s 1 b .
Consequently, we have shown that the functional V meq , s is nonnegative and vanishes if and only if θ s = θ s ^ and ρ s = ρ s ^ , which is in the case when the temperature and density fields reach their respective equilibrium values. The final formula for the solvent contribution reads
V meq , s = Ω ρ s c V , s θ s ^ θ s θ s ^ 1 ln θ s θ s ^ d v + Ω c V , s ( γ 1 ) θ s ^ ρ s ln ρ s ρ s ^ 1 b ρ s ^ 1 b ρ s ρ s ρ s ^ 1 b ρ s ^ d v .

6.4. Summary

We are now in a position to collect all partial results. The explicit formula for the functional (87) reads
V meq = Ω 1 2 ρ s v 2 d v + V meq , p + V meq , s = Ω 1 2 ρ s v 2 d v + k B θ s ^ Ω D M n p , θ s ^ φ M n p , θ s ^ ln φ M n p , θ s ^ φ M n p , θ s ^ + 1 d q d v + k B θ s ^ Ω n p ^ n p n p ^ ln n p n p ^ n p n p ^ + 1 d v + Ω ρ s c V , s θ s ^ θ s θ s ^ 1 ln θ s θ s ^ d v + Ω c V , s ( γ 1 ) θ s ^ ρ s ln ρ s ρ s ^ 1 b ρ s ^ 1 b ρ s ρ s ρ s ^ 1 b ρ s ^ d v
From the construction of the functional, it is clear that it decreases in time and that it is zero at equilibrium. In fact, we have an explicit formula for the time derivative of the functional:
d V meq d t = θ ^ Ω ξ d v .
The entropy production ξ can be identified from the right-hand side of (51), provided that we substitute for the fluxes j φ , x , j φ , q , j e and the Cauchy stress tensor T while using the formulae (52)–(54) and (62), which yields
ξ = 1 θ s 2 ν D δ : D δ + λ ˜ div x v 2 + 2 θ s ζ D 1 φ d U e d s s = 1 2 q q ref 2 + θ s θ ref d U η d s s = 1 2 q q ref 2 q q ref 2 φ + k B θ s q φ 2 d q + k B 2 θ s 2 ζ D 1 φ x φ 2 d q + κ x θ s 2 θ s 2 .
The first term in (117) monitors the approach of the velocity field to its equilibrium value (zero), the second term does the same regarding the configurational distribution function, the third term takes care of the polymer number density, the fourth term deals with the temperature field, and the last term monitors the approach of the density to the equilibrium density. However, we recall that we have not discussed the relation between the functional and a metric structure on the state space; hence, the functional does not, at the moment, deserve to be referred to as the Lyapunov functional. See also Dostalík et al. [50,51] for a similar application in the case of macroscopic viscoelastic rate-type models, where the authors have actually found a relation between the proposed Lyapunov type functionals and a metric on the state space.
Inspecting (119) and (118), we also note that the time derivative of V meq is proportional to the net entropy production that vanishes at the stationary spatially homogeneous equilibrium state (81). These observations conclude our rudimentary stability analysis.
Finally, we note that, on the right-hand side of (119), we explicitly see all entropy production mechanisms that are active in the given fluid. Because the entropy production terms, in fact, determine the fluxes in the sense of the discussion in Section 5.2, we see that the identification of the entropy production and the identification of the Helmholtz free energy indeed provide a complete characterisation of the fluid of interest.

7. Conclusions

We have provided a simple derivation of a model for non-isothermal flows of dilute polymeric fluids. The analysis outlined above shows that if the microscopic dynamics of the polymeric chains is governed by the Fokker–Planck equation with a centre-of-mass diffusion term, then it is possible to find a corresponding evolution equation for the temperature field in such a way that the whole system of governing equations is compatible with the first and second law of thermodynamics.
However, we note that, although the derivation outlined above seems to be completely self-consistent, some of the steps require insight into the microscopic dynamics of the polymeric chains. In particular, the interpretation of the coefficient ζ as hydrodynamic drag, see especially Formula (55) for the flux j φ , q , would be impossible without the insight into the microscopic nature of the model (this is however true also for the derivation proposed by Öttinger and Grmela [18]). Without such an insight, one could easily misidentify the flux j φ , q and lose the connection between the Fokker–Planck type equation and the underlying stochastic process at the level of individual polymeric chains.
While the presence of the centre-of-mass diffusion term is motivated by physical considerations, see for example El-Kareh and Leal [52], Bhave et al. [53], Schieber [54] and Degond and Liu [32], it turns out that its presence can be gainfully exploited in the mathematical analysis of the corresponding governing equations in the purely mechanical setting as well, see especially Barrett and Süli [31,33,55,56,57] or Feireisl et al. [34]. (See also Bulíček et al. [58], Bulíček et al. [59] and Bathory et al. [60] in the context of macroscopic models with a stress diffusion term.) One might hope that once thermodynamically consistent models for the full thermo-mechanical coupling are derived, similar rigorous mathematical results might also be obtained for the corresponding full system of governing equations, that is including the governing equation for the temperature.
In particular, the knowledge of the thermodynamic background leads to natural a priori estimates, and to the design of natural functionals that might be used to monitor the stability of steady equilibrium or non-equilibrium states or to analyse weak-strong uniqueness. (See for example Feireisl and Novotný [61], Feireisl et al. [62] and Feireisl and Pražák [63] for the discussion in the case of a compressible Navier–Stokes–Fourier fluid. In the context of purely mechanical models of dilute polymeric fluids and stability of the corresponding flows, see especially Jourdain et al. [64].) Indeed, we formally—that is under the assumption that the governing equations posses a strong solution—identify such a functional, which might be of use in such analyses.

Author Contributions

All four authors were equally involved in the formulation of the problem. The calculations were carried out by M.D. and V.P. The discussion of the results and the presentation of the material were carried out by all four authors. All authors have read and agreed to the published version of the manuscript.

Funding

J.M. thanks the Czech Science Foundation for the support through the project 18-12719S. M.D. and V.P. thank the Czech Science Foundation for the support through the project 20-11027X. M.D. has been also supported by Charles University Research program No. UNCE/SCI/023.

Acknowledgments

J.M. and V.P. acknowledge the membership of the Nečas Center for Mathematical Modeling (NCMM).

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Kramers, H.A. The behavior of macromolecules in inhomogeneous flow. J. Chem. Phys. 1946, 14, 415–424. [Google Scholar] [CrossRef]
  2. Bird, R.B.; Armstrong, R.C.; Hassager, O. Dynamics of Polymeric Liquids, 2nd ed.; Kinetic Theory; John Wiley & Sons: Hoboken, NJ, USA, 1987; Volume 2. [Google Scholar]
  3. Beris, A.N.; Edwards, B.J. Thermodynamics of Flowing Systems: With Internal Microstructure; Number 36 in Oxford Engineering Science Series; Oxford University Press: Oxford, UK, 1994. [Google Scholar]
  4. Öttinger, H.C. Stochastic Processes in Polymeric Fluids: Tools and Examples for Developing Simulation Algorithms; Springer: Berlin, Germay, 1996. [Google Scholar]
  5. Huilgol, R.R.; Phan-Thien, N. Fluid Mechanics of Viscoelasticity: General Principles, Constitutive Modelling, Analytical and Numerical Techniques; Rheology Series; Elsevier: Amsterdam, The Netherlands, 1997; Volume 6. [Google Scholar] [CrossRef]
  6. Dressler, M.; Edwards, B.J.; Öttinger, H.C. Macroscopic thermodynamics of flowing polymeric liquids. Rheol. Acta 1999, 38, 117–136. [Google Scholar] [CrossRef]
  7. Öttinger, H.C. Beyond Equilibrium Thermodynamics; John Wiley & Sons: Hoboken, NJ, USA, 2005. [Google Scholar]
  8. Kröger, M. Models for Polymeric and Anisotropic Liquids; Lecture Notes in Physics; Springer: Berlin, Germany, 2005; Volume 675. [Google Scholar]
  9. Lozinski, A.; Owens, R.G.; Phillips, T.N. The Langevin and Fokker–Planck equations in polymer rheology. In Numerical Methods for Non-Newtonian Fluids; Glowinski, R., Xu, J., Eds.; Handbook of Numerical Analysis; Elsevier: Amsterdam, The Netherlands, 2011; Volume 16, pp. 211–303. [Google Scholar] [CrossRef]
  10. Le Bris, C.; Lelièvre, T. Multiscale modelling of complex fluids: A mathematical initiation. In Multiscale Modeling and Simulation in Science; Engquist, B., Lötstedt, P., Runborg, O., Eds.; Lecture Notes in Computational Science and Engineering; Springer: Berlin/Heidelberg, Germany, 2009; Volume 66, pp. 49–137. [Google Scholar] [CrossRef] [Green Version]
  11. Lozinski, A.; Chauvière, C. A fast solver for Fokker-Planck equation applied to viscoelastic flows calculations: 2D FENE model. J. Comput. Phys. 2003, 189, 607–625. [Google Scholar] [CrossRef]
  12. Knezevic, D.J.; Süli, E. A heterogeneous alternating-direction method for a micro-macro dilute polymeric fluid model. M2AN Math. Model. Numer. Anal. 2009, 43, 1117–1156. [Google Scholar] [CrossRef] [Green Version]
  13. Mizerová, H.; She, B. A conservative scheme for the Fokker-Planck equation with applications to viscoelastic polymeric fluids. J. Comput. Phys. 2018, 374, 941–953. [Google Scholar] [CrossRef]
  14. Pavelka, M.; Klika, V.; Grmela, M. Multiscale Thermo-Dynamics; de Gruyter: Berlin, Germany, 2018. [Google Scholar]
  15. Peters, G.W.M.; Baaijens, F.P.T. Modelling of non-isothermal viscoelastic flows. J. Non-Newton. Fluid Mech. 1997, 68, 205–224. [Google Scholar] [CrossRef] [Green Version]
  16. Wapperom, P.; Hulsen, M.A. Thermodynamics of viscoelastic fluids: The temperature equation. J. Rheol. 1998, 42, 999–1019. [Google Scholar] [CrossRef] [Green Version]
  17. Hron, J.; Miloš, V.; Průša, V.; Souček, O.; Tůma, K. On thermodynamics of viscoelastic rate type fluids with temperature dependent material coefficients. Int. J. Non-Linear Mech. 2017, 95, 193–208. [Google Scholar] [CrossRef] [Green Version]
  18. Öttinger, H.C.; Grmela, M. Dynamics and thermodynamics of complex fluids. II. Illustrations of a general formalism. Phys. Rev. E 1997, 56, 6633–6655. [Google Scholar] [CrossRef]
  19. Coleman, B.D. On the stability of equilibrium states of general fluids. Arch. Ration. Mech. Anal. 1970, 36, 1–32. [Google Scholar] [CrossRef]
  20. Bulíček, M.; Málek, J.; Průša, V. Thermodynamics and stability of non-equilibrium steady states in open systems. Entropy 2019, 21, 704. [Google Scholar] [CrossRef] [Green Version]
  21. Málek, J.; Průša, V. Derivation of equations for continuum mechanics and thermodynamics of fluids. In Handbook of Mathematical Analysis in Mechanics of Viscous Fluids; Giga, Y., Novotný, A., Eds.; Springer: Berlin/Heidelberg, Germany, 2018; pp. 3–72. [Google Scholar] [CrossRef]
  22. Müller, I. Thermodynamics; Interaction of Mechanics and Mathematics; Pitman: London, UK, 1985. [Google Scholar]
  23. Gurtin, M.E.; Fried, E.; Anand, L. The Mechanics and Thermodynamics of Continua; Cambridge University Press: Cambridge, UK, 2010; p. xxii+694. [Google Scholar]
  24. Rajagopal, K.R.; Tao, L. Mechanics of Mixtures; Series on Advances in Mathematics for Applied Sciences; World Scientific Publishing Co. Inc.: River Edge, NJ, USA, 1995; Volume 35, p. xii+195. [Google Scholar]
  25. Hutter, K.; Jöhnk, K. Continuum Methods of Physical Modeling; Springer: Berlin, Germany, 2004. [Google Scholar]
  26. Souček, O.; Průša, V.; Málek, J.; Rajagopal, K.R. On the natural structure of thermodynamic potentials and fluxes in the theory of chemically non-reacting binary mixtures. Acta Mech. 2014, 225, 3157–3186. [Google Scholar] [CrossRef]
  27. Vinogradov, G.V.; Malkin, A.Y. Rheology of Polymers: Viscoelasticity and Flow of Polymers; Springer: Berlin/Heidelberg, Germany, 1980. [Google Scholar]
  28. Larson, R.G. Constitutive Equations for Polymer Melts and Solutions; Butterworths Series in Chemical Engineering; Butterworth-Heinemann: Oxford, UK, 1988. [Google Scholar] [CrossRef]
  29. Leonov, A.I.; Prokunin, A.N. Nonlinear Phenomena in Flows of Viscoelastic Polymer Fluids; Springer: Berlin/Heidelberg, Germany, 1994. [Google Scholar] [CrossRef]
  30. Tanner, R.; Walters, K. Rheology: An Historical Perspective; Rheology Series; Elsevier: Amsterdam, The Netherlands, 1998; Volume 7. [Google Scholar]
  31. Barrett, J.W.; Süli, E. Existence of global weak solutions to finitely extensible nonlinear bead-spring chain models for dilute polymers with variable density and viscosity. J. Diff. Equ. 2012, 253, 3610–3677. [Google Scholar] [CrossRef] [Green Version]
  32. Degond, P.; Liu, H. Kinetic models for polymers with inertial effects. Netw. Heterog. Media 2009, 4, 625–647. [Google Scholar] [CrossRef]
  33. Barrett, J.W.; Süli, E. Existence of global weak solutions to compressible isentropic finitely extensible bead-spring chain models for dilute polymers. Math. Model. Methods Appl. Sci. 2016, 26, 469–568. [Google Scholar] [CrossRef] [Green Version]
  34. Feireisl, E.; Lu, Y.; Süli, E. Dissipative weak solutions to compressible Navier–Stokes–Fokker–Planck systems with variable viscosity coefficients. J. Math. Anal. Appl. 2016, 443, 322–351. [Google Scholar] [CrossRef] [Green Version]
  35. Warner, H.R. Kinetic theory and rheology of dilute suspensions of finitely extendible dumbbells. Ind. Eng. Chem. Fundam. 1972, 11, 379–387. [Google Scholar] [CrossRef]
  36. Ericksen, J.L. Introduction to the Thermodynamics of Solids, revised ed.; Applied Mathematical Sciences; Springer: New York, NY, USA, 1998; Volume 131, p. xii+189. [Google Scholar] [CrossRef]
  37. Rajagopal, K.R.; Srinivasa, A.R. A thermodynamic frame work for rate type fluid models. J. Non-Newton. Fluid Mech. 2000, 88, 207–227. [Google Scholar] [CrossRef]
  38. Málek, J.; Rajagopal, K.R.; Tůma, K. On a variant of the Maxwell and Oldroyd-B models within the context of a thermodynamic basis. Int. J. Non-Linear Mech. 2015, 76, 42–47. [Google Scholar] [CrossRef]
  39. Málek, J.; Rajagopal, K.R.; Tůma, K. Derivation of the variants of the Burgers model using a thermodynamic approach and appealing to the concept of evolving natural configurations. Fluids 2018, 3, 69. [Google Scholar] [CrossRef] [Green Version]
  40. Málek, J.; Průša, V.; Skřivan, T.; Süli, E. Thermodynamics of viscoelastic rate-type fluids with stress diffusion. Phys. Fluids 2018, 30, 023101. [Google Scholar] [CrossRef]
  41. Dostalík, M.; Průša, V.; Skřivan, T. On diffusive variants of some classical viscoelastic rate-type models. AIP Conf. Proc. 2019, 2107, 020002. [Google Scholar] [CrossRef]
  42. Souček, O.; Orava, V.; Málek, J.; Bothe, D. A continuum model of heterogeneous catalysis: Thermodynamic framework for multicomponent bulk and surface phenomena coupled by sorption. Int. J. Eng. Sci. 2019, 138, 82–117. [Google Scholar] [CrossRef]
  43. Hatzikiriakos, S.G. Wall slip of molten polymers. Prog. Polym. Sci. 2012, 37, 624–643. [Google Scholar] [CrossRef]
  44. Coleman, B.D.; Greenberg, J.M. Thermodynamics and the stability of fluid motion. Arch. Ration. Mech. Anal. 1967, 25, 321–341. [Google Scholar] [CrossRef]
  45. Gurtin, M.E. Thermodynamics and the energy criterion for stability. Arch. Ration. Mech. Anal. 1973, 52, 93–103. [Google Scholar] [CrossRef]
  46. Gurtin, M.E. Thermodynamics and stability. Arch. Ration. Mech. Anal. 1975, 59, 63–96. [Google Scholar] [CrossRef]
  47. Grmela, M.; Öttinger, H.C. Dynamics and thermodynamics of complex fluids. I. Development of a general formalism. Phys. Rev. E 1997, 56, 6620–6632. [Google Scholar] [CrossRef]
  48. Le Métayer, O.; Saurel, R. The Noble–Abel stiffened-gas equation of state. Phys. Fluids 2016, 28, 046102. [Google Scholar] [CrossRef] [Green Version]
  49. Chiapolino, A.; Saurel, R. Extended Noble–Abel stiffened-gas equation of state for sub-and-supercritical liquid-gas systems far from the critical point. Fluids 2018, 3, 48. [Google Scholar] [CrossRef] [Green Version]
  50. Dostalík, M.; Průša, V.; Tůma, K. Finite amplitude stability of internal steady flows of the Giesekus viscoelastic rate-type fluid. Entropy 2019, 21, 1219. [Google Scholar] [CrossRef] [Green Version]
  51. Dostalík, M.; Průša, V.; Stein, J. Unconditional finite amplitude stability of a viscoelastic fluid in a mechanically isolated vessel with spatially non-uniform wall temperature. Math. Comput. Simul. 2020. [Google Scholar] [CrossRef]
  52. El-Kareh, A.W.; Leal, L.G. Existence of solutions for all Deborah numbers for a non-Newtonian model modified to include diffusion. J. Non-Newton. Fluid Mech. 1989, 33, 257–287. [Google Scholar] [CrossRef]
  53. Bhave, A.V.; Armstrong, R.C.; Brown, R.A. Kinetic theory and rheology of dilute, nonhomogeneous polymer solutions. J. Chem. Phys. 1991, 95, 2988–3000. [Google Scholar] [CrossRef]
  54. Schieber, J.D. Generalized Brownian configuration fields for Fokker–Planck equations including center-of-mass diffusion. J. Non-Newton. Fluid Mech. 2006, 135, 179–181. [Google Scholar] [CrossRef]
  55. Barrett, J.W.; Süli, E. Existence of global weak solutions to some regularized kinetic models for dilute polymers. Multiscale Model. Simul. 2007, 6, 506–546. [Google Scholar] [CrossRef] [Green Version]
  56. Barrett, J.W.; Süli, E. Existence and equilibration of global weak solutions to kinetic models for dilute polymers I: Finitely extensible nonlinear bead-spring chains. Math. Model. Methods Appl. Sci. 2011, 21, 1211–1289. [Google Scholar] [CrossRef] [Green Version]
  57. Barrett, J.W.; Süli, E. Existence and equilibration of global weak solutions to kinetic models for dilute polymers II: Hookean-type models. Math. Model. Methods Appl. Sci. 2012, 22, 1150024. [Google Scholar] [CrossRef] [Green Version]
  58. Bulíček, M.; Málek, J.; Průša, V.; Süli, E. A PDE-analysis for a class of thermodynamically compatible viscoelastic rate type fluids with stress diffusion. In Mathematical Analysis in Fluid Mechanics: Selected Recent Results; Danchin, R., Farwig, R., Neustupa, J., Penel, P., Eds.; Contemporary Mathematics; American Mathematical Society: Providence, RI, USA, 2018; Volume 710, pp. 25–53. [Google Scholar] [CrossRef] [Green Version]
  59. Bulíček, M.; Feireisl, E.; Málek, J. On a class of compressible viscoelastic rate-type fluids with stress-diffusion. Nonlinearity 2019, 32, 4665. [Google Scholar] [CrossRef] [Green Version]
  60. Bathory, M.; Bulíček, M.; Málek, J. Large Data Existence Theory for Three-Dimensional Unsteady Flows of Rate-Type Viscoelastic Fluids with Stress Diffusion. arXiv 2020, arXiv:2002.11224. [Google Scholar] [CrossRef]
  61. Feireisl, E.; Novotný, A. Weak–strong uniqueness property for the full Navier–Stokes–Fourier system. Arch. Ration. Mech. Anal. 2012, 204, 683–706. [Google Scholar] [CrossRef] [Green Version]
  62. Feireisl, E.; Jin, B.J.; Novotný, A. Relative entropies, suitable weak solutions, and weak-strong uniqueness for the compressible Navier-Stokes system. J. Math. Fluid Mech. 2012, 14, 717–730. [Google Scholar] [CrossRef]
  63. Feireisl, E.; Pražák, D. Asymptotic Behavior of Dynamical Systems in Fluid Mechanics; AIMS Series on Applied Mathematics; American Institute of Mathematical Sciences (AIMS): Springfield, MO, USA, 2010; Volume 4, p. xii+298. [Google Scholar]
  64. Jourdain, B.; Le Bris, C.; Lelièvre, T.; Otto, F. Long-time asymptotics of a multiscale model for polymeric fluid flows. Arch. Ration. Mech. Anal. 2006, 181, 97–148. [Google Scholar] [CrossRef]

Share and Cite

MDPI and ACS Style

Dostalík, M.; Málek, J.; Průša, V.; Süli, E. A Simple Construction of a Thermodynamically Consistent Mathematical Model for Non-Isothermal Flows of Dilute Compressible Polymeric Fluids. Fluids 2020, 5, 133. https://doi.org/10.3390/fluids5030133

AMA Style

Dostalík M, Málek J, Průša V, Süli E. A Simple Construction of a Thermodynamically Consistent Mathematical Model for Non-Isothermal Flows of Dilute Compressible Polymeric Fluids. Fluids. 2020; 5(3):133. https://doi.org/10.3390/fluids5030133

Chicago/Turabian Style

Dostalík, Mark, Josef Málek, Vít Průša, and Endre Süli. 2020. "A Simple Construction of a Thermodynamically Consistent Mathematical Model for Non-Isothermal Flows of Dilute Compressible Polymeric Fluids" Fluids 5, no. 3: 133. https://doi.org/10.3390/fluids5030133

APA Style

Dostalík, M., Málek, J., Průša, V., & Süli, E. (2020). A Simple Construction of a Thermodynamically Consistent Mathematical Model for Non-Isothermal Flows of Dilute Compressible Polymeric Fluids. Fluids, 5(3), 133. https://doi.org/10.3390/fluids5030133

Article Metrics

Back to TopTop