Kluitenberg-Verh\'as rheology of solids in the GENERIC framework

TL;DR

This paper reformulates the Kluitenberg-Verhás rheology model for solids within the GENERIC framework, enhancing theoretical understanding and potential practical applications in nonequilibrium thermodynamics.

Contribution

It presents the internal variable formulation of the Kluitenberg-Verhás model within the GENERIC framework, bridging rheology and thermodynamics.

Findings

Unified thermodynamic formulation of the Kluitenberg-Verhás model.

Enhanced understanding of solid rheology in nonequilibrium thermodynamics.

Potential for improved practical applications in material modeling.

Abstract

The internal variable methodology of nonequilibrium thermodynamics, with a symmetric tensorial internal variable, provides an important rheological model family for solids, the so-called Kluitenberg-Verh\'as model family [1]. This model family is distinguished not only from theoretical aspects but also on experimental grounds (see [2] for plastics and [3, 4, 5] for rocks). In this article, we present and discuss how the internal variable formulation of the Kluitenberg-Verh\'as model family can be presented in the nonequilibrium thermodynamical framework GENERIC (General Equation for the Non-Equilibrium Reversible-Irreversible Coupling) [6, 7, 8, 9], for the benefit of both thermodynamical methodologies as well as for promising practical applications.

Tables

Table 1. Table 1: How the four versions behave with respect to generic GENERIC expectations.

	variable $s$			variable $T$
	$l_A^d,l_A^s$ in $𝐌$	$l_A^d,l_A^s$ in $𝐋$		$l_A^d,l_A^s$ in $𝐌$	$l_A^d,l_A^s$ in $𝐋$
$𝐋$ fulfils Jacobi	✓	✓		✓	✓
$𝐌$ is symmetric	$\times$ if $l_A^d,l_A^s\ne 0$	✓		$\times$ if $l_A^d,l_A^s\ne 0$	✓

Equations

$$\dot{\mathbf{\epsilon }}={\left(\mathbf{v}\otimes \overleftarrow{\nabla }\right)}^{\mathrm{S}}.$$

$${\mathbf{x}}_0=\left(\mathbf{v},{\mathbf{\epsilon }}^{\mathrm{d}},{\mathbf{\epsilon }}^{\mathrm{s}},s_0\right).$$

$$\varrho \dot{\mathbf{v}}={\mathbf{\sigma }}_0\cdot \overleftarrow{\nabla }$$

$$\frac{\partial e_{\mathrm{int},0}}{\partial {\mathbf{\epsilon }}^{\mathrm{d}}}=\frac{1}{\varrho }{\mathbf{\sigma }}_0^{\mathrm{d}},\frac{\partial e_{\mathrm{int},0}}{\partial {\mathbf{\epsilon }}^{\mathrm{s}}}=\frac{1}{\varrho }{\mathbf{\sigma }}_0^{\mathrm{s}},\frac{\partial e_{\mathrm{int},0}}{\partial s_0}=T,$$

$$\mathbf{\sigma }$$

$$\varrho \dot{\mathbf{v}}$$

$$\displaystyle s={{s}_{0}}-\frac{1}{2}\operatorname{tr}{\left({\boldsymbol{\xi}}^{\rm d}{\boldsymbol{\xi}}^{\rm d}\right)}-\frac{1}{2}\operatorname{tr}{\Big{(}{\boldsymbol{\xi}}^{\rm s}{\boldsymbol{\xi}}^{\rm s}\Big{)}}.$$

$$\mathbf{x}=\left(\mathbf{v},{\mathbf{\epsilon }}^{\mathrm{d}},{\mathbf{\epsilon }}^{\mathrm{s}},s,{\mathbf{\xi }}^{\mathrm{d}},{\mathbf{\xi }}^{\mathrm{s}}\right),$$

$$\displaystyle{e}_{\rm int}\left(\mathbf{x}\right)={e}_{\rm int,0}\Big{(}{\boldsymbol{\varepsilon}}^{\rm d},{\boldsymbol{\varepsilon}}^{\rm s},{{s}_{0}}\big{(}s,{\boldsymbol{\xi}}^{\rm d},{\boldsymbol{\xi}}^{\rm s}\big{)}\Big{)}.$$

$$\varrho {\dot{e}}_{\mathrm{int}}=-{\mathbf{j}}_e\cdot \overleftarrow{\nabla }+\mathrm{tr}\left(\mathbf{\sigma }\dot{\mathbf{\epsilon }}\right),$$

$$\varrho {\dot{e}}_{\mathrm{int}}$$

$$\displaystyle=\operatorname{tr}\left({\boldsymbol{\sigma}}^{\rm d}_{0}{\dot{\boldsymbol{\varepsilon}}}^{\rm d}\right)+\operatorname{tr}\Big{(}{\boldsymbol{\sigma}}^{\rm s}_{0}{\dot{\boldsymbol{\varepsilon}}}^{\rm s}\Big{)}+\varrho T\dot{s}+\varrho T\operatorname{tr}\left({\boldsymbol{\xi}}^{\rm d}{\dot{\boldsymbol{\xi}}}\mskip 1.0mu^{\rm d}\right)+\varrho T\operatorname{tr}\Big{(}{\boldsymbol{\xi}}^{\rm s}{\dot{\boldsymbol{\xi}}}\mskip 1.0mu^{\rm s}\Big{)},$$

$$\displaystyle\varrho{\dot{e}}_{\rm int}=-\mathbf{j}_{e}\cdot{\overset{{\scriptscriptstyle\leftarrow}}{\nabla}}+\operatorname{tr}\left({\boldsymbol{\sigma}}^{\rm d}_{0}{\dot{\boldsymbol{\varepsilon}}}^{\rm d}\right)+\operatorname{tr}\Big{(}{\boldsymbol{\sigma}}^{\rm s}_{0}{\dot{\boldsymbol{\varepsilon}}}^{\rm s}\Big{)}+\operatorname{tr}\left({\hat{\boldsymbol{\sigma}}}^{\rm d}{\dot{\boldsymbol{\varepsilon}}}^{\rm d}\right)+\operatorname{tr}\Big{(}{\hat{\boldsymbol{\sigma}}}^{\rm s}{\dot{\boldsymbol{\varepsilon}}}^{\rm s}\Big{)}.$$

$$\displaystyle\varrho\dot{s}=-\frac{1}{T}\mathbf{j}_{e}\cdot{\overset{{\scriptscriptstyle\leftarrow}}{\nabla}}+\frac{1}{T}\operatorname{tr}\left({\hat{\boldsymbol{\sigma}}}^{\rm d}{\dot{\boldsymbol{\varepsilon}}}^{\rm d}\right)+\frac{1}{T}\operatorname{tr}\Big{(}{\hat{\boldsymbol{\sigma}}}^{\rm s}{\dot{\boldsymbol{\varepsilon}}}^{\rm s}\Big{)}-\varrho\operatorname{tr}\left({\boldsymbol{\xi}}^{\rm d}{\dot{\boldsymbol{\xi}}}\mskip 1.0mu^{\rm d}\right)-\varrho\operatorname{tr}\Big{(}{\boldsymbol{\xi}}^{\rm s}{\dot{\boldsymbol{\xi}}}\mskip 1.0mu^{\rm s}\Big{)}.$$

$$\varrho \dot{s}=-{\mathbf{j}}_s\cdot \overleftarrow{\nabla }+{\pi }_s$$

$$\displaystyle\pi_{s}=\varrho\dot{s}+\mathbf{j}_{s}\cdot{\overset{{\scriptscriptstyle\leftarrow}}{\nabla}}=\mathbf{j}_{e}\cdot\left(\frac{1}{T}\otimes{\overset{{\scriptscriptstyle\leftarrow}}{\nabla}}\right)+\frac{1}{T}\operatorname{tr}\Big{(}{\hat{\boldsymbol{\sigma}}}^{\rm d}{\dot{\boldsymbol{\varepsilon}}}^{\rm d}\Big{)}+\frac{1}{T}\operatorname{tr}\Big{(}{\hat{\boldsymbol{\sigma}}}^{\rm s}{\dot{\boldsymbol{\varepsilon}}}^{\rm s}\Big{)}-\varrho\operatorname{tr}\left({\boldsymbol{\xi}}^{\rm d}{\dot{\boldsymbol{\xi}}}\mskip 1.0mu^{\rm d}\right)-\varrho\operatorname{tr}\Big{(}{\boldsymbol{\xi}}^{\rm s}{\dot{\boldsymbol{\xi}}}\mskip 1.0mu^{\rm s}\Big{)}.$$

$${\hat{\mathbf{\sigma }}}^{\mathrm{d}}$$

$$\dot{\mathbf{\xi }}{}^{\mathrm{d}}$$

$$T{\pi }_s=\left(\begin{array}{cccc}{\dot{\mathbf{\epsilon }}}^{\mathrm{d}}& -\varrho T{\mathbf{\xi }}^{\mathrm{d}}& {\dot{\mathbf{\epsilon }}}^{\mathrm{s}}& -\varrho T{\mathbf{\xi }}^{\mathrm{s}}\end{array}\right)\left(\begin{array}{cccc}{l}_{11}^{\mathrm{d}}& l_{\mathrm{S}}^{\mathrm{d}}& 0& 0\\ l_{\mathrm{S}}^{\mathrm{d}}& l_{22}^{\mathrm{d}}& 0& 0\\ 0& 0& l_{11}^{\mathrm{s}}& l_{\mathrm{S}}^{\mathrm{s}}\\ 0& 0& l_{\mathrm{S}}^{\mathrm{s}}& l_{22}^{\mathrm{s}}\end{array}\right)\left(\begin{array}{c}{\dot{\mathbf{\epsilon }}}^{\mathrm{d}}\\ -\varrho T{\mathbf{\xi }}^{\mathrm{d}}\\ {\dot{\mathbf{\epsilon }}}^{\mathrm{s}}\\ -\varrho T{\mathbf{\xi }}^{\mathrm{s}}\end{array}\right)$$

$$l_{11}^{\mathrm{d}}$$

$$l_{11}^{\mathrm{s}}$$

$${\hat{\mathbf{\sigma }}}^{\mathrm{d}}$$

$$\dot{\mathbf{\xi }}{}^{\mathrm{d}}$$

$${\mathbf{\sigma }}^{\mathrm{d}}+{\tau }^{\mathrm{d}}{\dot{\mathbf{\sigma }}}^{\mathrm{d}}$$

$${\mathbf{\sigma }}^{\mathrm{s}}+{\tau }^{\mathrm{s}}{\dot{\mathbf{\sigma }}}^{\mathrm{s}}$$

$$\displaystyle\frac{\mathrm{d}\mathbf{x}}{\mathrm{d}t}=\mathbf{L}\mathopen{}\mathclose{{}\left(\mathbf{x}}\right)\frac{\delta E}{\delta\mathbf{x}}+\mathbf{M}\mathopen{}\mathclose{{}\left(\mathbf{x}}\right)\frac{\delta S}{\delta\mathbf{x}},$$

$$\mathbf{L}=-{\mathbf{L}}^{\mathrm{T}}$$

$$\mathbf{M}\frac{\delta E}{\delta \mathbf{x}}=0,$$

$$\mathbf{L}\frac{\delta S}{\delta \mathbf{x}}=0$$

$$\left\{A,\left\{B,C\right\}\right\}+\left\{B,\left\{C,A\right\}\right\}+\left\{C,\left\{A,B\right\}\right\}=0$$

$$\left\{A,B\right\}:={\int }_V\frac{\delta A}{\delta \mathbf{x}}\mathbf{L}\frac{\delta B}{\delta \mathbf{x}}\mathrm{d}V.$$

$$\mathbf{M}={\mathbf{M}}^{\mathrm{T}},$$

$$\left[A,B\right]:={\int }_V\frac{\delta A}{\delta \mathbf{x}}\mathbf{M}\frac{\delta B}{\delta \mathbf{x}}\mathrm{d}V$$

Full text

1 Introduction

The internal variable approach of nonequilibrium thermodynamics, with a symmetric tensorial internal variable, provides a distinguished model family – the Kluitenberg–Verhás model family [1] (covering the Hooke, Kelvin–Voigt, Maxwell, Poynting–Thomson and Jeffrey models as special cases) – for the rheology of solids. This family is significant not only from theoretical perspective but also for experimental applications [2, 3, 4, 5]. GENERIC (General Equation for Nonequilibrium Reversible–Irreversible Coupling) is an attractive general framework for nonequilibrium thermodynamical models (see, e.g., [6, 7, 8, 9]). Whenever a new nonequilibrium thermodynamical model emerges, it is advantageous and recommended to check how it suits the frame of GENERIC. Here, we investigate how the internal variable formulation leading to the Kluitenberg–Verhás model family can be represented in GENERIC.

For the main part of the paper, specific entropy is treated as one of the state variables—a choice natural from principal aspects. Later, in an alternative version, temperature is used, instead—which formulation may be more convenient for certain engineering applications.

We believe that the relationship between the internal variable framework and GENERIC may be fruitful for both approaches, providing

• •

insight concerning the theoretical side,

• •

wider applicability,

• •

suggestions for novel numerical methods (see, e.g., [10, 11] for such a promising direction), and

• •

beneficial connection of such numerical approaches with analytical results (e.g., [12]).

In this respect, this paper intends to serve as a case study.

Notably, certain aspects of our treatment could be presented at some more general, systematic and methodological level (for example, performing the deviatoric–spherical separation of tensors in a multiplicative split form [13]). Here, we follow a simple direct approach.

2 Necessary elements I: The internal variable formulation of

rheology of solids

We start with a summary and generalization of the internal variable approach to the Kluitenberg–Verhás rheological model family of solids [1]. The discussion is generalized in that the derivation in [1] neglected thermal expansion and started from Hookean elasticity, while the version here is free of those restrictions, only isotropy of the material being assumed.

The small-strain regime is considered, where strain $\boldsymbol{\varepsilon}$ is small ($\lvert{\boldsymbol{\varepsilon}}\rvert\ll 1$), there is no need to distinguish spacetime and material manifold variables and vectors/tensors—accordingly, aspects of objectivity and spacetime compatibility [14, 15] are not addressed here—, mass density $\varrho$ is constant, and one can relate time derivative (partial coinciding with substantial) of strain with the symmetrized gradient111${\overset{{\scriptscriptstyle\leftarrow}}{\nabla}}$ and ${\overset{{\scriptscriptstyle\rightarrow}}{\nabla}}$ act to the left and to the right, respectively, reflecting proper tensorial order. of the velocity field $\mathbf{v}$,

[TABLE]

Due to the isotropy of the material, the deviatoric–spherical decomposition of symmetric tensors plays here an important role (the spherical part of, *e.g., *strain is proportional to the identity tensor, ${\boldsymbol{\varepsilon}}^{\rm s}=\frac{1}{3}\left(\operatorname{tr}\boldsymbol{\varepsilon}\right)\mathbf{1}$, while ${\boldsymbol{\varepsilon}}^{\rm d}=\boldsymbol{\varepsilon}-{\boldsymbol{\varepsilon}}^{\rm s}$ is its deviatoric part). With ${{s}_{0}}$ denoting mass-specific entropy, our variables will be

[TABLE]

The balance of linear momentum is

[TABLE]

with the divergence of stress ${\boldsymbol{\sigma}}_{0}$ on the rhs, where ${\boldsymbol{\sigma}}_{0}$ is related to a (mass-)specific internal energy ${e}_{\rm int,0}\left({\boldsymbol{\varepsilon}}^{\rm d},{\boldsymbol{\varepsilon}}^{\rm s},{{s}_{0}}\right)$ as shown by the partial derivatives222$\frac{\partial{e}_{\rm int,0}}{\partial{\boldsymbol{\varepsilon}}^{\rm d}}$ is purely deviatoric since, for any two tensors $\mathbf{A},\mathbf{B}$, $\,{{A}^{\rm s}}_{ij}{{B}^{\rm d}}_{ij}=0\,$ so \,\frac{\partial{e}_{\rm int,0}}{\partial{\varepsilon}^{\rm d}_{ij}}\mathrm{d}{\varepsilon}^{\rm d}_{ij}={\Big{(}\frac{\partial{e}_{\rm int,0}}{\partial{\varepsilon}^{\rm d}_{ij}}\Big{)}\mskip-1.0mu}^{\rm d}\mathrm{d}{\varepsilon}^{\rm d}_{ij}\,.\,

[TABLE]

$T$ standing for temperature.

Rheology is a behaviour most manifest in the mechanical aspect so, for an internal variable description of it, in conform with that strain and stress are second order symmetric tensors, we introduce a symmetric tensorial internal variable $\boldsymbol{\xi}$. Mechanical effects of rheology are to be embodied by a $\boldsymbol{\xi}$ dependent extension of stress:

[TABLE]

We conceive rheology as irreversibility-related so specific entropy is also assumed to be influenced; concavity concerns combined with Morse’s lemma for smooth enough $\boldsymbol{\xi}$ dependence and nonzero second derivative in $\boldsymbol{\xi}$ conclude in the variable transformation ${{s}_{0}}\to s$ [16]

[TABLE]

Correspondingly, specific internal energy expressed in terms of the extended collection of variables,

[TABLE]

is of the form

[TABLE]

The balance of internal energy is

[TABLE]

where $\mathbf{j}_{e}$ denotes heat current density, and the only source term considered is related to mechanical power. Substituting (4), (5) and (7) into (10), on the one hand we obtain

[TABLE]

and on the other hand

[TABLE]

The rhs of (11) is to be equal to the rhs of (12), which leads to

[TABLE]

Since the balance of the extended entropy is to be of the form

[TABLE]

with entropy current density $\mathbf{j}_{s}$ chosen to be the usual $\,\mathbf{j}_{s}=\frac{1}{T}\mathbf{j}_{e}\,$, and entropy production $\pi_{s}$, in the light of (13), we can write

[TABLE]

Positive semidefiniteness of entropy production can be ensured for the first term via $\,\,\mathbf{j}_{e}=\lambda\left(\frac{1}{T}\otimes{\overset{{\scriptscriptstyle\leftarrow}}{\nabla}}\right)\,,\,\,$ $\lambda\geq 0$ (Fourier heat conduction, a vectorial part that cannot isotropically couple to the remaining, tensorial, terms; hence, for simplicity, heat conduction is neglected in what follows), and via Onsagerian equations concerning the further terms, with independent deviatoric and spherical parts because of isotropy:

[TABLE]

with appropriate conditions on the deviatoric coefficients $l_{ij}^{{\rm d}}$ and the spherical ones $l_{ij}^{{\rm s}}$, each of which are going to be assumed constant for simplicity333No principal difficulties are expected when these coefficients are $\mathbf{x}$ dependent, and here we intend to keep formulae relatively short.. These conditions can be read off from the quadratic form obtained by substituting (16)–(17) into (15), which yields444The upper right and lower left two-by-two submatrices contain only zero elements due to the isotropic decoupling of deviatoric and spherical parts.

[TABLE]

with $l^{\rm d}_{\rm S}=\frac{1}{2}\left(l_{12}^{{\rm d}}+l_{21}^{{\rm d}}\right)$ and $l^{\rm s}_{\rm S}=\frac{1}{2}\left(l_{12}^{{\rm s}}+l_{21}^{{\rm s}}\right)$. Hence, the four-by-four coefficient matrix in the middle is required to be positive semidefinite, which necessitates for the coefficients, using Sylvester’s criteria,

[TABLE]

We remark that, both in (19) and (20), the three conditions are not independent: the third one and either of the first two ones imply the remaining one. It is important to notice that the antisymmetric part of the coefficient matrix does not contribute to entropy production. We can emphasize this by dividing the Onsagerian equations (16)–(17) into two parts:

[TABLE]

with \,l^{\rm d}_{\rm A}=\frac{1}{2}\mathopen{}\mathclose{{}\left(l_{12}^{{\rm d}}-l_{21}^{{\rm d}}}\right)\, and \,l^{\rm s}_{\rm A}=\frac{1}{2}\Big{(}l_{12}^{{\rm s}}-l_{21}^{{\rm s}}\Big{)}\,.

It is to be noted that, in general, the coefficient matrices $l^{{\rm d}}$, $l^{{\rm s}}$ need not be symmetric nor antisymmetric, corresponding to that the concrete physical interpretation of $\boldsymbol{\xi}$ may not be available and the behaviour of $\boldsymbol{\xi}$ under time reflection might not be purely sign preserving/flipping.555Rocks are one example of complex enough materials that may require such a description.

To see that this model family covers classic rheological models like Kelvin–Voigt and Poynting–Thomson, one can start with the special case of Hooke elasticity, and eliminating the internal variable leads, in the isothermal approximation (constant $l_{11}^{\rm d,s}$, $l_{\rm A}^{\rm d,s}$, $\varrho Tl_{\rm S}^{\rm d,s}$, $\varrho Tl_{22}^{\rm d,s}$), to the Kluitenberg–Verhás model family [1],

[TABLE]

with necessary and sufficient thermodynamical inequality conditions on the coefficients $\tau^{\rm d,s},E_{0}^{\rm d,s},E_{1}^{\rm d,s},E_{2}^{\rm d,s}$ stemming from (19)–(20) (for further details on the elimination and the inequalities, see [1], Section 2.3).

3 Necessary elements II: Summary of the GENERIC framework

In GENERIC [6, 7], time evolution of the collection of state variables (fields, in case of continuum models like ours here), $\mathbf{x}$, is formulated as

[TABLE]

where the operator matrix $\mathbf{L}$ acts on the column vector that is the functional derivative of the energy functional $E$ of $\mathbf{x}$, and the operator matrix $\mathbf{M}$ acts on the column vector that is the functional derivative of the entropy functional $S$ of $\mathbf{x}$. $\mathbf{L}$ is required to be antisymmetric,

[TABLE]

(T denoting transpose which, for operators, means not merely matrix index transposition but includes operator adjoint). Thanks to this and the degeneracy condition

[TABLE]

energy is conserved, $\,\frac{\mathrm{d}E}{\mathrm{d}t}=0\,$. In parallel, the other degeneracy requirement

[TABLE]

ensures that the first term on the rhs of (25) does not increase entropy, and $\mathbf{M}$ is demanded to be positive semidefinite to lead to $\,\frac{\mathrm{d}S}{\mathrm{d}t}\geq 0\,$ eventually. That $\mathbf{L}\frac{\delta E}{\delta\mathbf{x}}$ is related to reversible dynamics is manifested further by also prescribing the Jacobi identity

[TABLE]

($A$, $B$, $C$ arbitrary functionals) for the bilinear generalized Poisson bracket

[TABLE]

Consequently, the first term on the rhs of (25) can be interpreted as a reversible—generalized Hamiltonian—time evolution contribution (vector field) while the second term (another vector field) embodies irreversible time evolution contribution to dynamics.

Analogously to (26) and (30), imposing symmetricity for $\mathbf{M}$,

[TABLE]

induces that the bilinear product

[TABLE]

is positive semidefinite, $\left[A,A\right]\geq 0$. The latter bracket (32) completes the former one (30) in the sense that time evolution for any functional $A$ can be expressed as

[TABLE]

A constructive and productive way to generate the irreversible contribution to dynamics is to derive it from a dissipation potential [7, 17, 18]. Assuming a dissipation potential is, on the other side, not necessary and reduces the level of generality of the GENERIC framework [19].

4 Internal variable rheology of solids realized in the GENERIC

formulation

Section 2 has actually been given in a form to provide preparations for the present one, where we establish GENERIC form for the $\boldsymbol{\xi}$-described rheology of solids. The set of variables $\mathbf{x}$ is (8), the energy functional consists of the internal energy contribution (9) supplemented by the kinetic energy related one, and the entropy functional is straightforward:

[TABLE]

The corresponding functional derivatives are

[TABLE]

The nontrivial task is to identify $\mathbf{L}$ and $\mathbf{M}$. Concerning the time evolution of the state variables, we know (6), (1), (13) and (16)–(17) so we conjecture the decomposition to reversible and irreversible parts as666Beware that, if we prefer to write $\mathbf{x}$ and $\dot{\mathbf{x}}$ as column vectors then (36) should contain row vectors, as being covectors with respect to the vector space of $\mathbf{x}$. However, then $\mathbf{L}$ and $\mathbf{M}$ could not be displayed as customary square matrices. The misleading double meaning of column vectors could be resolved by writing $\mathbf{x}$-covectors like in (36) as column vectors but within $[\>\>]$ instead of $(\>\>)$. Then, correspondingly, $\mathbf{L}$, $\mathbf{M}$ were to be written within $(\>\>]$. Here, we decided not to use this convention but at least to draw attention to that distinction between vectors and covectors is not only principally important but also avoids considerable confusion during calculations.

[TABLE]

The governing principle for this decision for decomposition is that, since dissipation and irreversibility are related to entropy production and to the internal variable, the reversible vector field should not contain them but only pure fluid mechanics.

Then $\mathbf{L}$ can directly read off from the first term on the rhs of (37):

[TABLE]

with $\bullet$ denoting the ‘slot’ where the operator acts. This $\mathbf{L}$ apparently fulfils the degeneracy condition $\,\mathbf{L}\frac{\delta S}{\delta\mathbf{x}}=\mathbf{0}\,$.

To prove antisymmetry of $\mathbf{L}$, let us take the corresponding bracket (30):

[TABLE]

where $A$ and $B$ are arbitrary functionals of the state variables, and abbreviations of the kind

[TABLE]

have been introduced. Using indices (with Einstein convention and the Kronecker symbol $\delta_{ij}$), we have

[TABLE]

leading to

[TABLE]

Via integration by parts and omitting surface terms,

[TABLE]

Now, taking into consideration that the functional derivative of a scalar functional with respect to a symmetric tensor is symmetric, with respect to a deviatoric tensor is deviatoric, and with respect to a spherical tensor is spherical, the first term in the integrand can be reformulated as

[TABLE]

and the second term can be treated analogously. Next, it is easy to show that the third term is \,B_{\mathbf{v}}\cdot\mathopen{}\mathclose{{}\left(A_{{\boldsymbol{\varepsilon}}^{\rm d}}\cdot{\overset{{\scriptscriptstyle\leftarrow}}{\nabla}}}\right)\,. Further, in the terms that contain $\delta_{ij}$, the $j$ index can be changed to $i$, hence, in these terms we find the gradient of the trace of a tensor. Therefore, the fourth term contains trace of a deviatoric tensor—which is traceless—so this term gives zero contribution. Finally, the fifth term contains a spherical tensor and thus can be rewritten as \,\frac{1}{3}\partial_{k}\mathopen{}\mathclose{{}\left(A_{{\boldsymbol{\varepsilon}}^{\rm s}}}\right)_{ii}\mathopen{}\mathclose{{}\left(B_{\mathbf{v}}}\right)_{k}=\partial_{k}\mathopen{}\mathclose{{}\left(A_{{\boldsymbol{\varepsilon}}^{\rm s}}}\right)_{ik}\mathopen{}\mathclose{{}\left(B_{\mathbf{v}}}\right)_{i}\,. To sum up, we find

[TABLE]

antisymmetry is revealed.

Since $\mathbf{L}$ is independent of the state variables and is antisymmetric, the generalized Poisson bracket also satisfies the Jacobi identity [6].

Now let us turn towards the irreversible side: the operator matrix $\mathbf{M}$ can be constructed from the second term of the time evolution equation (37) and the degeneracy condition (27); we find

[TABLE]

with777An outline and order of obtaining the components is as follows: $M_{24}=0$ and $M_{34}=0$ from (37), $M_{42}=0$ and $M_{43}=0$ from (31), $M_{56}=0$ and $M_{65}=0$ from (37), $M_{44}$ from (37), $M_{54}$ and $M_{64}$ from (37), $M_{45}$ and $M_{46}$ from (31), $M_{41}$ from (27), $M_{14}$ from (31), $M_{51}$ and $M_{55}$ from (37), $M_{61}$ and $M_{66}$ from (37), $M_{15}$ and $M_{16}$ from (31), $M_{11}$ from (27), the unconstrained remaining components can be set to zero.

[TABLE]

One can notice that this $\mathbf{M}$ is not symmetric—see the elements that contain $l_{12}^{{\rm d}}$ vs. $l_{21}^{{\rm d}}$, as well as the ones with $l_{12}^{{\rm s}}$ vs. $l_{21}^{{\rm s}}$. As mentioned in Section 2, when we have no additional microscopic or experimental information about $\boldsymbol{\xi}$ and about the corresponding coefficients $l^{{\rm d}}$, $l^{{\rm s}}$ then we cannot exclude that antisymmetric parts $l^{\rm d}_{\rm A}$, $l^{\rm s}_{\rm A}$ appear in the dynamics.

On the other side, positive semidefiniteness can be shown by reformulating the integrand of the irreversible bracket $\left[A,A\right]$ to a quadratic expression. More closely, we can form a matrix that contains nonnegative elements and the Onsagerian coefficients, and (19)–(20) just prove to be the conditions that ensure positive semidefiniteness. The calculation is straightforward but lengthy.

Actually, the whole realization of $\boldsymbol{\xi}$-based rheology provided above is straightforward (if lengthy), and is expected to work for nonconstant coefficient matrices $l^{{\rm d}}$, $l^{{\rm s}}$ as well. However, specifically for constant coefficients, an alternative version is also possible: implementing the antisymmetric part of the coefficient matrices, that is, the constants $l^{\rm d}_{\rm A}$, $l^{\rm s}_{\rm A}$, in the reversible part of the time evolution.

For this case, let us use the prepared (21)–(22) form of the Onsagerian equations. Rearranging the time evolution equation is straightforward, and we find for the alternative reversible operator matrix $\mathbf{L}^{\prime}$

[TABLE]

while the elements of the alternative irreversible operator matrix $\mathbf{M}^{\prime}$ are very similar to (48): we just have to change all $l_{12}^{{\rm d}}$ and $l_{21}^{{\rm d}}$ to $l^{\rm d}_{\rm S}$ and, similarly, $l_{12}^{{\rm s}}$ and $l_{21}^{{\rm s}}$ to $l^{\rm s}_{\rm S}$.

The degeneracy criteria, antisymmetry of $\mathbf{L}^{\prime}$, the Jacobi identity for the generalized Poisson brackets, and positive semidefiniteness for $\mathbf{M}^{\prime}$ prove to be satisfied. Moreover, in this case the symmetric property of $\mathbf{M}^{\prime}$ also holds.

We repeat that this alternative realization is valid only for constant Onsagerian coefficients as otherwise the Jacobi identity were violated.

This latter variant appears rather counter-intuitive from principal point of view but may be beneficial for numerical solutions, e.g., to have as much symplectic part in the numerical scheme as possible—see [10, 11] about the importance of this.

5 Temperature as state variable

For mechanical engineering applications and evaluations of experiments (see, e.g., [2]), it can be beneficial to use temperature, instead of entropy, as state variable. Then the collection of state variables is

[TABLE]

To keep the discussion as concrete and basic as possible, let us choose the simplest constitutive equation for the internal energy, linear in temperature with constant specific heat $c$, containing elastic energy related to Hooke’s law,

[TABLE]

and neglecting thermal expansion—which is manifested in the separation of strain dependence and temperature dependence—:

[TABLE]

Temperature has the same relationship to specific entropy as seen in Section 4, now utilized in the opposite direction (*i.e., *what is a variable and what is a function): The thermodynamical consistency condition $\,\frac{\mathrm{d}{{\tilde{s}}_{0}}}{\mathrm{d}T}=\frac{1}{T}\frac{\mathrm{d}{\tilde{e}}_{\rm th}}{\mathrm{d}T}\,$ that follows from the Gibbs relation [and which is the manifestation of the third equation of (4)] leads to

[TABLE]

with auxiliary constants $\tilde{s}_{{\rm aux}}$, $T_{{\rm aux}}$, and the extension (7) induces

[TABLE]

or, expressing temperature,

[TABLE]

Now the energy and entropy functionals are

[TABLE]

with corresponding functional derivatives

[TABLE]

We perform a transformation of variables $\,\mathbf{x}\to\tilde{\mathbf{x}}\,$, to which the transformation (operator) matrix $\,\mathbf{Q}=\frac{\delta\tilde{\mathbf{x}}}{\delta\mathbf{x}}\,$ is accompanied. This $\mathbf{Q}$ can be used to establish the relationship between the original and transformed reversible and irreversible operator matrices [6]:

[TABLE]

In the present current special case, we change only the fourth state variable (from $s$ to $T$), so only the fourth row of $\mathbf{Q}$ contains nontrivial elements. Furthermore, since (55) does not contain nonlocal (gradient) terms, we can realize the transformation directly in the form

[TABLE]

where $\mathbf{I}$ denotes the fourth order identity tensor [the identity that maps tensors to tensors, *i.e., *themselves]. Then, using (59) yields888Results (61), (62) can also be derived directly from the time evolution formula and the degeneracy conditions.

[TABLE]

so all the requirements of GENERIC—antisymmetry, Jacobi identity, degeneracy—prove to hold for $\tilde{\mathbf{L}}$, and we find that the structure of $\tilde{\mathbf{M}}$ is the same as of $\mathbf{M}$ [see (47)], with the elements

[TABLE]

The variable transformation is expected to preserve the structure of GENERIC ([6], page 26 in Section 1.2.4). Indeed, by substituting (52) into (10) and rewriting it in terms of temperature, the evolution equation for $T$ is obtained, and turns out to coincide with the fourth row of $\tilde{\mathbf{M}}\frac{\delta\tilde{S}}{\delta\tilde{\mathbf{x}}}$ so the whole evolution equation has been preserved under the transformation. Meanwhile, as in the variable $s$ case $\mathbf{M}$ has proved nonsymmetric for nonzero $l^{\rm d}_{\rm A}$ or $l^{\rm s}_{\rm A}$, $\tilde{\mathbf{M}}$ behaves the same way.

Now let us repeat moving the $l^{\rm d}_{\rm A}$ and $l^{\rm s}_{\rm A}$ related part of the dynamics to the reversible part. We find the antisymmetric

[TABLE]

with

[TABLE]

and that $\tilde{\mathbf{M}}^{\prime}$ can be obtained from $\tilde{\mathbf{M}}$ like $\mathbf{M}^{\prime}$ from $\mathbf{M}$, *i.e., *changing all $l_{12}^{{\rm d}}$ and $l_{21}^{{\rm d}}$ to $l^{\rm d}_{\rm S}$ and $l_{12}^{{\rm s}}$ and $l_{21}^{{\rm s}}$ to $l^{\rm s}_{\rm S}$. Symmetricity and positive semidefiniteness of $\tilde{\mathbf{M}}$, the degeneracy criteria, as well as the Jacobi identity related to $\tilde{\mathbf{L}}^{\prime}$ are all satisfied.

It is to be noted here that, while the Jacobi property of $\tilde{\mathbf{L}}^{\prime}$ is foreseen on general grounds – any variable transformation is expected to preserve the Jacobi identity ([6], page 26 in Section 1.2.4) –, checking it directly is nonstraightforward. The difficulty is related to the task of identifying total divergences of multiple products among the numerous terms. When we used the application jacobi.m [20] – with appropriately increased memory limit and run-time limit –, it could not confirm the Jacobi identity (while it found its validity for $\mathbf{L}$, $\mathbf{L}^{\prime}$ and $\tilde{\mathbf{L}}$ seamlessly). Instead, we have verified the Jacobi property of $\tilde{\mathbf{L}}^{\prime}$ both by hand and via an own symbolic code. A key element was an advantageous convention for classifying and grouping terms, which has reduced the number of terms from thousands to hundreds, enabling to observe the remaining cancellations.

6 Conclusions

The results can be summarized according to Table 1.

As a task for the future, a finite deformation version would be welcome. How deeply this will require to address objectivity and spacetime compatibility aspects [14, 15] is an open question.

In parallel, the current small-strain version could be numerically (*e.g., *along the lines of [10, 11]) applied for concrete problems. For example, the recent analytical results in [12] promise methodologically useful outcome since those considerations are done in the force equilibrial approximation [zero lhs in (6), an approximation valuable for various engineering situations], which is a challenge for GENERIC with its explicite time evolution formulation. Principal as well as numerically working solutions to this compelling situation can provide fruitful contribution to both science and engineering.

Acknowledgement

The authors thank Markus Hütter for his assistance concerning the initial steps, Hans Christian Öttinger, Miroslav Grmela, and Michal Pavelka for valuable remarks, the participants of the workshop IWNET 2018 for the useful discussions, and Karl Heinz Hoffmann and the Reviewers for helpful suggestions.

The work was supported by the Hungarian grants National Research, Development and Innovation Office – NKFIH 116375 and NKFIH KH 130378, and by the Higher Education Excellence Program of the Ministry of Human Capacities in the frame of Nanotechnology research area of Budapest University of Technology and Economics (BME FIKP-NANO).