1. Introduction
Granular materials are assemblies of solid particles that interact through inelastic collisions and contacts. The thermodynamic responses of granular materials have wide-ranging applications due to their ubiquity in industry and nature, where the vast majority consist of non-spherical particles. In industrial settings the efficient, accurate and safe transportation, sorting and storing of granular materials are required. Also, all naturally occurring granular flows (landslides, avalanches and bedload sediment) are composed of non-spherical particles (sand particles, rocks or pebbles) where the prediction, prevention and control of such flows are critical. Beyond natural and industrial relevance, granular physics present a fascinating scientific challenge since it can exhibit complex responses that can span from solids to gases that are governed by the microscopic interactions. Therefore, the development of reliable granular thermodynamic models is essential.
Extensive research has been devoted to granular thermodynamics, ranging from gas phase (Haff Reference Haff1983; Mehta Reference Mehta2007), fluid phase (MiDi Reference MiDi2004; Jop, Forterre & Pouliquen Reference Jop, Forterre and Pouliquen2006; Sun & Sundaresan Reference Sun and Sundaresan2011) to quasistatic solid phase (Alonso-Marroquin & Herrmann Reference Alonso-Marroquin and Herrmann2004; García-Rojo et al. Reference García-Rojo, Alonso-Marroquín and Herrmann2005). However, most existing studies neglect the additional complexities introduced by the shape and orientation of non-spherical particles. For instance, the kinetic theory of granular gases, where the interaction between particles is built around binary collisions between hard spherical particles (Jenkins & Savage Reference Jenkins and Savage1983; Jenkins & Berzi Reference Jenkins and Berzi2010), does not account for the distinct responses arising from the shape- and orientation-dependent collisions between particles.
Figure 1. A granular gas of oriented particles depicting the density, velocity, thermodynamic temperature, particle agitations and orientations.
Depicted in figure 1 is a granular gas consisting of identical oriented particles and their relevant physical quantities. The figure highlights these quantities, including density (showing regions of low and high density), velocity (arrows indicating velocity magnitude and direction) and thermodynamic temperature (showing cold and hot particles), which are standard quantities not unique to granular gases. For granular gases, the particle agitation including both translational and rotational motions (showing low and high particle agitation) is an additional physical quantity. The intensity of the agitations is represented by a scalar quantity named granular temperature or velocity fluctuations. The associated kinetic energy is an additional form of internal energy named granular energy or kinetic energy fluctuations. In addition, for granular systems of oriented particles, the orientation quantity represented by the orientation tensor, is the statistical distribution of particle orientations, ranging from isotropic to aligned configurations, as shown in figure 1.
The physical behaviour of granular systems with densities below the critical density, i.e. no long-term contacts between particles are established, is governed solely by particle interactions through collisions, and such systems are therefore referred to as granular gases. The frequency of the collisions is proportional to the granular agitation, that is, the granular temperature. The macroscopic behaviour of granular gases is directly related to the microscopic interactions between the particles, which involve multiple aspects. During a collision, particles exchange momentum, thermal energy and granular energy. For inelastic collisions, part of the granular energy is also transferred into thermal energy. For oriented particles, these interactions are affected by the shape and orientations of the particles. Consequently, collisions also lead to an exchange of orientation between particles.
In granular gases, the particles’ thermal energy is typically not explicitly considered. While the dissipation, which is the rate of transfer of granular energy into thermal energy due to inelastic interactions, is explicitly considered in the balance of granular energy, the balance of thermal energy is typically ignored since the duration of collisions is too short to exchange thermal energy. However, this exchange of thermal energy cannot be neglected for dense granular systems where long-term contacts are established.
In the following, we are interested in the macroscopic behaviour of granular gases rather than the microscopic discrete behaviour. We adopt a continuum approach where the discrete properties of the particles are replaced by continuum fields. However, the continuum description is motivated by the underlying interactions between particles.
It is well known that the state of particle agitations, quantified by the granular temperature, substantially influences the macroscopic response of granular gases and it is an essential component in formulating the governing equations of such granular systems (Jenkins & Savage Reference Jenkins and Savage1983; Jenkins & Richman Reference Jenkins and Richman1985; Garzó & Dufty Reference Garzó and Dufty1999; Goldhirsch Reference Goldhirsch2003, Reference Goldhirsch2008; Jenkins & Berzi Reference Jenkins and Berzi2012). Its significant role in granular systems of oriented particles is also reported in Berzi et al. (Reference Berzi, Thai-Quang, Guo and Curtis2016, Reference Berzi, Thai-Quang, Guo and Curtis2017) and Pol, Artoni & Richard (Reference Pol, Artoni and Richard2023).
Here, we consider a system of particles with identical shape and size. Oriented particles are particles with sufficient symmetry, such as cylinders, disks, rice, lentils and similar shapes. A particle’s orientation is identified by an axis vector
$\pm k_i$
, where the difference between
$+ k_i$
and
$- k_i$
is not detectable. Then, an oriented particle is defined as a particle whose orientation can be fully characterized by the dyad
$k_ik_j$
, thereby removing the arbitrariness associated with the sign of
$\pm k_i$
. Particles such as cones and pyramids are not sufficiently symmetric, since there is a detectable difference between
$+k_i$
and
$-k_i$
, making their orientation more complex, hence, the dyad
$k_ik_j$
may not be sufficient to represent their orientation.
It is recognized that the shape and orientation of particles play a significant role in the response of granular systems consisting of oriented particles. It is observed that orientation is strongly affected by the flow which is coupled to the rheological properties that are governed by the orientation (Berzi et al. Reference Berzi, Thai-Quang, Guo and Curtis2016; Nagy et al. Reference Nagy, Claudin, Börzsönyi and Somfai2017; Nadler, Guillard & Einav Reference Nadler, Guillard and Einav2018; Trulsson Reference Trulsson2018; Pol et al. Reference Pol, Artoni and Richard2023) using discrete numerical simulations. These observations are also supported by experiments using advanced experimental techniques that provide important insights into the flow, orientation and rheology. Image processing and X-ray tomography were applied to measure more complex flows, such as split-bottom, silo, hopper, down incline, bottlenecks and confined cylindrical domain (Stannarius et al. Reference Stannarius, Wegner, Szabó and Börzsönyi2013; Börzsönyi et al. Reference Börzsönyi, Somfai, Szabó, Wegner, Mier, Rose and Stannarius2016; Guillard, Marks & Einav Reference Guillard, Marks and Einav2017; Parisi et al. Reference Parisi, Cruz Hidalgo and Zuriguel2018; Hidalgo et al. Reference Hidalgo, Szabó, Gillemot, Börzsönyi and Weinhart2018; Pol et al. Reference Pol, Artoni, Richard, da, Nunes and Gabrieli2022) clearly identifying a strong coupling between the orientation and the flow. Studying inhomogeneous problems (Hidalgo et al. Reference Hidalgo, Szabó, Gillemot, Börzsönyi and Weinhart2018; Pol et al. Reference Pol, Artoni, Richard, da, Nunes and Gabrieli2022; Berzi, Vescovi & Nadler Reference Berzi, Vescovi and Nadler2025) revels the properties of orientation flux and orientation boundary conditions (Amereh & Nadler Reference Amereh and Nadler2022; Berzi et al. Reference Berzi, Vescovi and Nadler2025).
There are many similarities between liquid crystals (Gennes & Prost Reference Gennes and Prost1993) and oriented granular materials, particularly the treatment to the directionality. However, the underlying physics of the two systems is fundamentally different. In granular materials, interactions through collisions play a central role, while electromagnetic interactions that are essential in liquid crystals are absent.
In (Leslie Reference Leslie1968) liquid crystal are taken as micropolar fluids (Ericksen Reference Ericksen1961; Leslie Reference Leslie1966) where a director is identified with the preferred direction of the molecules. They can be treated as anisotropic fluids (Hand Reference Hand1962) with a symmetric second order-structural tensor that can represent a more general orientation than a single director. The structural state variable
$A_{\textit{ij}}$
used here is equivalent to the Q-tensor appearing in the free energy density of liquid crystals (De Gennes Reference De Gennes1969; Mottram & Newton Reference Mottram and Newton2014).
We adopt the view of Truesdell & Noll (Reference Truesdell, Noll and Fl’ugge1965): ‘Of course, physical theory must be based on experience, but experiment comes after, rather than before, theory. Without theoretical concepts one would neither know what experiments to perform nor be able to interpret their outcome’. This statement refers to phenomenological continuum approach which is the direction we pursue in this work.
Therefore, before we consider the available data from discrete simulations and experiments, we must first establish the theoretical framework for granular gases of oriented particles. To this end, we use the principles of non-equilibrium thermodynamics (De Groot & Mazur Reference De Groot and Mazur1962) to rigorously formulate a general framework for such systems. Only then, the framework is invigorated by determining the transport coefficients of granular systems of oriented particles using the available data. It should be noted that even the granular kinetic theory (Jenkins & Savage Reference Jenkins and Savage1983; Jenkins & Richman Reference Jenkins and Richman1985; Garzó & Dufty Reference Garzó and Dufty1999; Soto Reference Soto2004; Jenkins & Berzi Reference Jenkins and Berzi2012), which provides very important results for spherical particles, must be phenomenologically adjusted to obtain good agreement with discrete simulations and experiments of realistic granular systems. Furthermore, for oriented particles due to the additional complexity, the phenomenological approach cannot be avoided.
Non-equilibrium thermodynamics (De Groot & Mazur Reference De Groot and Mazur1962) is a continuum approach that systematically formulate the theory in terms of transport equations, in balance law form, and phenomenological transport coefficients. The first step is to choose the state variables which are the relevant continuum quantities describing the state of the system. For simple gas, these are typically the mass density, velocity and thermal temperature (Struchtrup Reference Struchtrup2024). For granular gases the granular temperature is an additional state variable (Garzó & Dufty Reference Garzó and Dufty1999; Jenkins & Berzi Reference Jenkins and Berzi2012). For granular gases consisting of identical spherical particles, the mass density, velocity, thermal temperature and granular temperature are sufficient to describe the state of the system. For granular gases consisting of identical oriented particles, the orientation tensor (Nadler et al. Reference Nadler, Guillard and Einav2018) is an additional state variable.
The conservation laws for mass and linear momentum describe the evolution for the density and velocity fields (De Groot & Mazur Reference De Groot and Mazur1962). The balance of internal energy is decomposed into the balance of thermal energy and the balance of granular energy, where energy exchange occurs through inelastic interaction between particles. For oriented particles, the balance of the orientation tensor is postulated. The conservation and balance laws include non-convective fluxes and productions for which constitutive relations must be determined to close the system of equations.
The challenging task of formulating these constitutive laws is done systematically in non-equilibrium thermodynamics (De Groot & Mazur Reference De Groot and Mazur1962) by considering the balance of entropy and the second law of thermodynamics. The formulation of the fluxes and productions must satisfy the second law of thermodynamics as well as the objectivity requirement (Truesdell & Noll Reference Truesdell, Noll and Fl’ugge1965). A careful consideration of the balance of entropy and the Gibbs equation (Müller Reference Müller1985) provides the expression for the entropy flux and production. For granular gases of oriented particles, the Gibbs equation is extended to include the contributions of the granular energy and orientation to entropy. The entropy production, which appears as the product of thermodynamical forces and fluxes, allows us to identify the thermodynamical forces, i.e. the deviations from equilibrium, and the associated fluxes that drive the system towards equilibrium. Since the second law of thermodynamics requires that the entropy production must be non-negative and vanish at equilibrium, phenomenological relationships for the fluxes in terms of the thermodynamical forces are established.
The formulation becomes systematic using classical linear non-equilibrium thermodynamics (De Groot & Mazur Reference De Groot and Mazur1962; Struchtrup Reference Struchtrup2024) where the fluxes are considered to be linear in the thermodynamic forces. The Onsager–Casimir reciprocal relations (Onsager Reference Onsager1931a , Reference Onsagerb ; Casimir Reference Casimir1945) provide further restrictions that simplify the formulation. Finally, constitutive laws for the fluxes using the representation theorems (Wang Reference Wang1969) that satisfy the second law of thermodynamics are formulated in terms of the phenomenological transport coefficients that can be a function of the state variables.
Only now when the theoretical framework, which is derived from basic physical principles, is established, it is ready to correlate with available results and experimental data from physical experiments or discrete microscopic simulations. This approach provides guidelines to what is required, allowed and, even more important, what is not allowed, when formulating constitutive laws for granular gases of oriented particles. Next, the phenomenological transport coefficients are determined to obtain satisfying agreement with results from the granular kinetic theory and experimental data.
Finally, the closed set of transport equations can be applied to solve initial-boundary valued problems typically using numerical methods to solve the resulting system of partial deferential equations. However, solving relevant initial-boundary value problems is left to future work as it is not the focus of this work.
The outline of the paper is as follows: § 2 defines the state variables and formulates the associated conservation laws, balance laws and the second law of thermodynamics. The constitutive laws are presented in § 3. In § 4 the initial and boundary conditions are formulated. Applications to granular systems of oriented particles and comparison with results available in the literature are presented in § 5. Section 6 is a summary of the governing equations. Finally, conclusions are presented in § 7.
Throughout the paper, we use index notation to represent tensors. Scalars (tensors of order zero) have no indices, vectors (tensors of order zero) have one index, second-order tensors have two indices, etc. Symmetric pairs of indices are denoted by parentheses,
$()$
, antisymmetric pairs of indices by brackets
$[]$
, and symmetric and traceless pairs of indices by angle brackets
$\langle \rangle$
, so that, e.g.
\begin{align} \varPsi _{i(jklm)n} &= \frac {1}{2}\left ( \varPsi _{ijklmn} + \varPsi _{imkljn} \right )\! , \end{align}
\begin{align} \varPsi _{i[jklm]n} &= \frac {1}{2}\left ( \varPsi _{ijklmn} - \varPsi _{imkljn} \right )\! , \end{align}
\begin{align} \varPsi _{i\langle jklm\rangle n} &= \frac {1}{2}\left ( \varPsi _{ijklmn} + \varPsi _{imkljn} \right ) - \frac {1}{3}\varPsi _{ipklpn}{\delta }_{jm} , \end{align}
where the summation convention, a summation over the range of an index that appears twice in a term, is implied unless otherwise indicated. Also, nested brackets are evaluated from the inside outwards. The concise notation for powers of second-order tensors
\begin{equation} A^2_{\textit{ij}} = A_{il}A_{\textit{lj}} , A^3_{\textit{ij}} = A_{il}A_{\textit{lm}}A_{\textit{mj}} \end{equation}
is used. The common notation
\begin{equation} \varPsi _{i,j} = \frac {{\partial } \varPsi _{i}}{{\partial } x_{i}} \end{equation}
for spatial partial derivatives is used.
2. State variables, conservation law, balance laws and the second law of thermodynamics
It is sufficient for a system consisting of a single isotropic material to take the density
$\rho$
, velocity
$v_i$
and thermodynamic temperature
$\theta$
as the state variables to fully describe the state of the system (De Groot & Mazur Reference De Groot and Mazur1962; Struchtrup Reference Struchtrup2024). For granular systems of oriented particles, we must extend the state variables to include the granular temperature
$T$
and the orientation tensor
$A_{\textit{ij}}$
, both of which are essential to characterize the equilibrium and non-equilibrium states (Vescovi, Nadler & Berzi Reference Vescovi, Nadler and Berzi2024).
The evolutions of the state variables
\begin{equation} \{\rho , v_i, \theta , T, A_{\textit{ij}}\} \end{equation}
are governed by the conservation of mass and momentum, and the balance laws for the thermal energy, granular energy and the orientation, which are discussed next.
2.1. Balance laws of oriented granular gases
A balance law is a partial differential equation that has the general form
\begin{equation} \rho \dot {\varGamma } + \gamma _{l,l} = \digamma , \end{equation}
where
$\varGamma$
is a state variable,
$\dot {\varGamma }={\partial }\varGamma /{\partial } t + v_l\varGamma ,_l$
denotes the material time derivative,
$\gamma _l$
is the non-convective flux of
$\varGamma$
and
$\digamma$
is the production. If the production
$\digamma$
vanishes identically for all processes it is a conservation law rather than a balance law.
2.1.1. Conservation of mass and linear momentum
The governing equations of the evolution for the density,
$\rho$
, and the velocity,
$v_i$
, are the standard (De Groot & Mazur Reference De Groot and Mazur1962; Struchtrup Reference Struchtrup2024) conservation laws for mass and linear momentum
\begin{align} \dot {\rho } + \rho v_{l,l} &=0 , \end{align}
\begin{align} \rho \dot {v}_i + t_{il,l} &= 0 , \end{align}
where
$t_{\textit{ij}}$
is the symmetric pressure tensor, satisfying the balance of angular momentum.
2.1.2. Balance of internal energy
The specific internal energy,
$E=u+U$
, is decomposed to thermal energy,
$u$
, and granular energy,
$U$
. The thermal energy is the internal energy of the solid particles taking the standard form
\begin{equation} u=c_v \theta , \end{equation}
where
$c_v$
is the specific heat. The granular energy,
$U$
, represents the kinetic energy associated with the particle velocity fluctuations
\begin{equation} U = \frac {5}{2} T , \end{equation}
where the coefficient (5/2), as opposed to the classical (3/2), used for spherical particles, captures the contributions of the five active degrees of freedom (three translational and two rotational) using the equipartition principle that energy is shared equally between all degrees of freedom, similarly to diatomic gases (Nagnibeda & Kustova Reference Nagnibeda and Kustova2009). The rotation about the particle’s symmetry axis is considered negligible, as it is activated primarily by surface friction rather than by collisions.
We take the granular temperature as a scalar quantity that characterizes the overall magnitude of the granular energy that includes the translational and rotational fluctuations. This is justified by the argument that, for oriented particles, both types of fluctuations contribute to interparticle collisions and energy exchange. This contrasts with spherical particles, for which only translational velocity fluctuations lead to collisions.
Also, the solid particles are taken to be rigid, and thus deformation energy is neglected. The thermodynamic temperature,
$\theta$
, (hereafter referred to as `temperature’), and the granular temperature,
$T$
, are expressed in energy unit for convenience, i.e. with dimensions of velocity squared,
$\rm(length/time)^2$
, hence the specific heat,
$c_v$
, is dimensionless.
In non-equilibrium states, thermal and granular energies are transported through their respective fluxes, driving the system towards equilibrium. Additionally, energy exchange occurs between granular and thermal modes due to inelastic interactions between particles, as well as through internal mechanical power.
The balance of thermal energy and granular energy are
\begin{align} \rho \dot {u} + q_{l,l} &= r , \end{align}
\begin{align} \rho \dot {U} + Q_{l,l} &= -r-t_{\textit{lm}} v_{(l,m)} , \end{align}
where
$q_{i}$
is the thermal heat flux,
$r$
denotes the rate of energy exchange from granular to thermal energy due to inelastic particle interactions,
$Q_i$
is the granular heat flux and
$t_{\textit{lm}} v_{(l,m)}$
is the internal mechanical power.
2.1.3. Balance of orientation
The statistical microscopic distribution of particle orientations (Hand Reference Hand1962; De Gennes Reference De Gennes1969; Mottram & Newton Reference Mottram and Newton2014; Nadler et al. Reference Nadler, Guillard and Einav2018) is represented by a second-order tensor
\begin{equation} A_{\textit{ij}} = \oint _{s^2} f(\boldsymbol{k})k_i k_j {\rm d}a - \frac {1}{3}\delta _{\textit{ij}} , \end{equation}
where
$A_{\textit{ij}}$
is a symmetric traceless second-order tensor,
$s^2$
denotes a unit sphere,
$f(\boldsymbol{k})$
is the probability density function in the direction
$\boldsymbol{k}$
with
$\oint _{s^2} f(\boldsymbol{k}){\rm d}a=1$
, and
$k_i$
are the components of the axis vector
$\boldsymbol{k}$
. The orientation tensor represents the deviation from a uniform (isotropic) distribution,
$f(\boldsymbol{k})=(4\pi )^{-1} \mapsto A_{\textit{ij}}=0$
, which is taken to be the equilibrium state. The limit state away from the equilibrium occurs when all particles align in a single direction, say
$\boldsymbol{k}^0$
, such that
$f(\boldsymbol{k}) = {\delta }(\boldsymbol{k}-\boldsymbol{k}^0) \mapsto A_{\textit{ij}}=k^0_i k^0_j-1/3 {\delta }_{\textit{ij}}$
, in which case the eigenvalues are
$\{2/3,-1/3,-1/3\}$
. It follows that the eigenvalues of
$A_{\textit{ij}}$
are bounded by
$[-1/3,2/3]$
. The linear invariant of the orientation tensor vanishes identically,
$A_{ll}=0$
, and the two nonlinear invariants (Wang Reference Wang1969) can be defined as
\begin{equation} i_1 = 3/2 A^2_{ll} , \quad i_2 = 9/2 A^3_{ll} , \end{equation}
which are scaled such that they vanish at isotropic distribution and attain the maximum value of one under perfect alignment with
$i_1 \in [0,1]$
and
$i_2 \in [-1/8,1]$
.
The balance of orientation is given by
\begin{equation} \rho{A}_{\textit{ij}} + G_{\textit{ijl},l} = P_{\textit{ij}} , \end{equation}
where
${A}_{\textit{ij}}=\dot {A}_{\textit{ij}} + A_{il} v_{[l,j]} - v_{[i,l]} A_{\textit{lj}}$
is the objective Jaumann–Zaremba material time derivative (Noll Reference Noll1968),
$G_{\textit{ijk}}$
is the non-convective orientation flux and
$P_{\textit{ij}}$
is the orientation production. The orientation flux
$G_{\textit{ijk}}$
and the orientation production
$P_{\textit{ij}}$
are constrained to ensure that the eigenvalues of the orientation tensor remain within their bounds.
These constraints can be formulated in terms of the eigenvalues using the spectral representation of the orientation tensor
\begin{equation} A_{\textit{ij}} = \sum _{l=1}^3 A^l u^l_i u^l_j , \end{equation}
where
$\{A^l\}$
are the eigenvalues and
$\{u^l_i\}$
are the associated orthonormal eigenvectors such that
$u^i_lu^j_l={\delta }_{\textit{ij}}$
. The eigenvalues have the properties
$\sum _{l=1}^3 A^l =0$
and
$A^l \in [-1/3,2/3]$
. The evolution of an eigenvalue is then expressed by
\begin{equation} \rho \dot {A}^i = \left ( P_{\textit{lm}} - G_{lmn,n} \right ) u^i_l u^i_m , \end{equation}
with no summation over the superscript
$i$
. The restriction requires that when an eigenvalue reaches its minimum value
$-1/3$
, its time derivative is non-negative, or equivalently, when an eigenvalue reaches its maximum value
$2/3$
, its time derivative is non-positive. This can be formulated as
\begin{equation} \left[P_{\textit{lm}}u^i_l u^i_m - G_{lmn,n}u^i_l u^i_m\right]_{A^i=-1/3} \geqslant 0 . \end{equation}
Since the orientation flux and production are independent, each term must satisfy the restriction, individually
\begin{align} \left[P_{\textit{lm}}u^i_l u^i_m\right]_{A^i=-1/3} &\geqslant 0 , \end{align}
\begin{align} \left[G_{lmn,n}u^i_l u^i_m\right]_{A^i=-1/3} &\leqslant 0 . \end{align}
These constraints must be considered when constructing the constitutive laws from the orientation flux,
$G_{\textit{ijk}}$
and production
$P_{\textit{ij}}$
.
Equations (2.3)–(2.11) form the governing system of equations for the evolution of the state variables. To close the system, constitutive laws for the fluxes and productions must be developed. For this, we apply the principles of linear irreversible thermodynamics.
2.2. Entropy of oriented granular gas
The collisions between particles leads to exchange of momentum and energy resulting a collective behaviour of the granular system which is described by few macroscopic quantities. Such a collective behaviour is entropic as it is driven by processes that maximize the entropy.
The fluxes, which depend on the system and the non-equilibrium state, must be provided for closure of the balance laws. We use the term `flux’ broadly to refer to any quantity that drives the system towards equilibrium including the standard fluxes
$\{q_i, Q_i, t_{\textit{ij}}, G_{\textit{ijk}}\}$
and productions
$\{r, P_{\textit{ij}}\}$
. The second law of thermodynamics states that entropy can only be created, ensuring directionality of irreversible processes. We introduce the balance of entropy
\begin{equation} \rho \dot {s} + \phi _{l,l} = \varSigma \ \geqslant 0 , \end{equation}
where
$s$
is the specific entropy,
$\phi _i$
is the entropy flux and
$\varSigma$
is the non-negative entropy production.
The processes considered here are such that local interaction between particles through collisions are sufficiently frequent that local equilibrium is obtained at much shorter time scales than global equilibrium. It follows that it is required that collisions are weakly inelastic, that is, small deviations from elastic collisions. Hence, we assume a separation of time scales such that processes are in local equilibrium but in global non-equilibrium. The global non-equilibrium state must be maintained by external excitation.
Entropy is an additive quantity, such that at local equilibrium the differential form of the Gibbs equation is universal in the sense that it holds without specifying the equations of state. The dependence of entropy on the variables is given through the Gibbs equation (Struchtrup Reference Struchtrup2024), which we write for granular gases by extending the standard Gibbs equation to include the additional granular contributions, granular energy and orientation, to the entropy as
\begin{equation} {\rm d}s = - p {\rm d}\rho /\bigl(T \rho ^2\bigr) + {\rm d}u/\theta + {\rm d}U/T - cA_{\textit{lm}} \,{\rm d}A_{\textit{lm}} . \end{equation}
Here,
$p$
is the pressure and the last two terms are the additional granular contributions. The last term accounts for the dependence of the entropy on the granular orientation, which we propose to be quadratic in the orientational tensor
$A_{\textit{ij}}$
similarly to the first term equation (2) in De Gennes (Reference De Gennes1969). Indeed, entropy (2.18) must be concave and maximum at equilibrium, where
$A_{\textit{ij}}=0$
, which required that
$c\geqslant 0$
. Due to the chosen temperature scale, the entropy is dimensionless, hence the parameter
$c$
is dimensionless as well, but it could depend on the state variables. We emphasize that such state dependence is restricted by the requirement that the Gibbs equation (2.18) be an exact differential, i.e. an entropy function exists. In the present work, we consider constitutive choices for which these integrability conditions are satisfied, a specific case of which is taking
$c$
to be constant.
Substitution of (2.3)–(2.11) into the Gibbs equation (2.18) yields
\begin{equation} \rho \dot {s} = \bigl(r-q_{l,l}\bigr)/\theta - \left(r+Q_{l,l} - t^v_{\textit{lm}} v_{(l,m)}\right)/T - cA_{\textit{lm}} \bigl(P_{\textit{lm}}-G_{lmn,n} \bigr) , \end{equation}
where the viscous stress is defined as
\begin{equation} t^{v}_{\textit{ij}} = p {\delta }_{\textit{ij}} - t_{\textit{ij}} . \end{equation}
Constructing the divergence terms
\begin{align} q_{l,l}/\theta &= (q_{l}/\theta )_{,l} - q_{l} (1/\theta )_{,l} , \end{align}
\begin{align} Q_{l,l}/T &= (Q_{l}/T)_{,l} - Q_{l} (1/T)_{,l} , \end{align}
\begin{align} A_{\textit{lm}} G_{lmn,n} &= (G_{lmn} A_{\textit{lm}})_{,n} - G_{lmn} A_{lm,n} , \end{align}
and substitution into (2.19) results in the entropy balance law of the form (2.17)
\begin{align} \rho \dot {s} + \left (q_{i}/\theta + Q_{i}/T - c A_{\textit{lm}} G_{lmi}\right )_{,i} &=( 1/\theta - 1/T ) r + q_{l} (1/\theta )_{,l} + Q_{l} (1/ T)_{,l} \nonumber \\ &\quad + t^v_{\textit{lm}} v_{(l,m)}/T - c A_{\textit{lm}} P_{\textit{lm}} -c A_{lm,n} G_{lmn}, \end{align}
from which we can identify the entropy flux and production.
By (2.24) the entropy flux is
\begin{equation} \phi _i = q_{i}/\theta + Q_{i}/T - c A_{\textit{lm}} G_{lmi} , \end{equation}
with contributions related to thermal energy, granular energy and orientation fluxes. And the entropy production is
\begin{equation} \varSigma = ( 1/\theta - 1/T ) r + q_{l} (1/\theta )_{,l} + Q_{l} (1/ T)_{,l} + t^v_{\textit{lm}} v_{(l,m)}/T - c A_{\textit{lm}} P_{\textit{lm}} -c A_{lm,n} G_{lmn} , \end{equation}
which includes the standard entropy productions due to thermal heat flux and viscous irreversibility (De Groot & Mazur Reference De Groot and Mazur1962; Struchtrup Reference Struchtrup2024), as well as non-standard contributions arising from granular-–thermal energy exchange, granular heat flux, orientation production and orientation flux.
The second law of thermodynamics requires non-negative entropy production,
$\varSigma \geqslant 0$
, which is used in the next section to construct constitutive laws for the fluxes
$\{r,q_i,Q_i,t^v_{\textit{ij}},P_{\textit{ij}},G_{\textit{ijk}}\}$
.
3. Constitutive laws
Equation (2.26) expresses the entropy production as the sum of products of forces and fluxes. The forces
$\mathcal{F}_\alpha$
represent deviations from equilibrium and the fluxes
$\mathcal{J}_\alpha$
drive the system towards equilibrium, a state of maximum entropy. Hence, the non-negative entropy production imposes the requirement that
\begin{align} \varSigma =\sum _{\alpha } \mathcal{F}_\alpha \mathcal{J}_\alpha \geqslant 0 \, . \end{align}
From (2.26) the forces and fluxes are identified as
\begin{align} \mathcal{F}_\alpha &= \begin{bmatrix} \; 1/\theta - 1/T , ( 1/\theta )_{,i} , ( 1/T )_{,i} , v_{(i,j)} , - c A_{\textit{ij}} , - c A_{\textit{ij},k} \; \end{bmatrix} \!, \end{align}
\begin{align} \mathcal{J}_\alpha &= \begin{bmatrix} \; r , q_i , Q_i , t^v_{\textit{ij}}/T , P_{\textit{ij}} , G_{\textit{ijk}} \; \end{bmatrix} \!, \end{align}
respectively. At equilibrium, all forces
$\mathcal{F}_\alpha$
vanish, hence the associated fluxes
$\mathcal{J}_\alpha$
must also vanish to preserve the equilibrium state. We follow the procedures of classical linear non-equilibrium thermodynamics (De Groot & Mazur Reference De Groot and Mazur1962; Struchtrup Reference Struchtrup2024) where fluxes are linear functions of the forces with coefficients that depend on the local state. These linear relations ensure preservation of equilibrium and can be constraint to satisfy non-negative entropy production for all processes.
The linear relations are written as
\begin{align} \mathcal{J}_\alpha =\sum _{\beta } L_{\alpha \beta } \mathcal{F}_\beta , \end{align}
where
$L_{\alpha \beta }$
is the matrix of phenomenological coefficients representing the constitutive properties of the system. Since both forces and fluxes are Galilean-invariant (objective), the matrix of phenomenological coefficients must also be Galilean-invariant, thus can only depend on objective state variables
\begin{equation} \{\rho , \theta , T , A_{\textit{ij}}\} , \end{equation}
hence, excluding,
$v_i$
, the velocity. These coefficients must be determined from microscopic theories or from data obtained from microscopic simulations or physical experiments.
Substituting (3.4) into (3.1), the entropy production reads
\begin{align} \varSigma = \sum _{\alpha } \mathcal{F}_\alpha \mathcal{J}_\alpha = \sum _{\alpha } \sum _{\beta } \mathcal{F}_\alpha L_{\alpha \beta } \mathcal{F}_\beta \geqslant 0 \,, \end{align}
which implies that
$L_{\alpha \beta }$
must be a positive semidefinite matrix.
It should be noted that the Curie principle, stating that fluxes and forces of different tensorial orders do not couple (De Groot & Mazur Reference De Groot and Mazur1962), only holds for isotropic systems, and does not hold here due to the anisotropy property,
$A_{\textit{ij}}$
, which introduces coupling between tensors of different orders. Furthermore, the unusual dual role of
$A_{\textit{ij}}$
as both a state variable and a force introduces additional complexity.
Considering all possible couplings between tenors of different orders, the general form of the phenomenological matrix is
\begin{equation} \left [\begin{array}{c} r \\[3pt] q_{i} \\[3pt] Q_{i} \\[3pt] t^{v}_{\textit{ij}}/T \\[3pt] P_{\textit{ij}}\\[3pt] G_{\textit{ijk}} \end{array}\right ]=\left [\begin{array}{cccccc} B^{1} & 0 & 0 & B_{\textit{lm}}^{2} & B_{\textit{lm}}^{3} & 0\\[3pt] 0 & C_{il}^{1} & C_{il}^{2} & 0 & 0 & C_{ilmn}^{3}\\[3pt] 0 & D_{il}^{1} & D_{il}^{2} & 0 & 0 & D_{ilmn}^{3}\\[3pt] E_{\textit{ij}}^{1} & 0 & 0 & E_{\textit{ijlm}}^{2} & E_{\textit{ijlm}}^{3} & E_{ijlmn}^{4}\\[3pt] F_{\textit{ij}}^{1} & 0 & 0 & F_{\textit{ijlm}}^{2} & F_{\textit{ijlm}}^{3} & F_{ijlmn}^{4}\\[3pt] 0 & H_{\textit{ijkl}}^{1} & H_{\textit{ijkl}}^{2} & H_{\textit{ijklm}}^{3} & H_{\textit{ijklm}}^{4} & H_{ijklmn}^{5} \end{array}\right ]\left [\begin{array}{c} 1/\theta -1/ T \\[3pt] (1/\theta )_{,l} \\[3pt] (1/T)_{,l} \\[3pt] v_{(l,m)} \\[3pt] - c A_{\textit{lm}} \\[3pt] - c A_{lm,n} \end{array}\right ] , \end{equation}
where the phenomenological coefficients are functions of the state variables (3.5). The vanishing terms are due to the unavailability of a vector and a third-order tensor with the required symmetry that could link the corresponding terms.
Additional simplification can be obtained by utilizing the Onsager–Casimir reciprocal relations (Onsager Reference Onsager1931a , Reference Onsagerb ; Casimir Reference Casimir1945) stating that macroscopic equations are invariant under time reversal. Although the microscopic dynamics of granular gases are not time-reversible due to inelastic collisions, the particles considered here correspond to small deviations from elastic collisions. In this limit, it has been shown (Garzó Reference Garzó2019) that reciprocal relations hold exactly for elastic collisions and remain a good approximation for weakly inelastic systems. On this basis, their use is justified as an approximation, not as a fundamental requirement. In what follows, the reciprocal relations are adopted as an approximation that significantly simplifies the structure of the phenomenological transport coefficients.
The reciprocity relationship (De Groot & Mazur Reference De Groot and Mazur1962) implies that
$L_{\alpha \beta }=L_{\beta \alpha }$
, if
$\mathcal{F}_\alpha$
and
$\mathcal{F}_\beta$
exhibit the same transformation under time-reversal, otherwise
$L_{\alpha \beta }=-L_{\beta \alpha }$
. Specifically, by consideration of (3.2), only
$\mathcal{F}_4=v_{(i,j)}$
transforms under time-reversal, hence,
$L_{\alpha \beta }=L_{\beta \alpha }$
except
$L_{\alpha 4}=-L_{4\alpha }$
. Applying the reciprocal relations to (3.7) yields the final from of the phenomenological coefficient matrix
\begin{equation} \left [\begin{array}{c} r\\ q_{i}\\[3pt] Q_{i}\\[3pt] t_{\textit{ij}}^{v}/T\\[3pt] P_{\textit{ij}}\\[3pt] G_{\textit{ijk}} \end{array}\right ]=\left [\begin{array}{cccccc} \mathcal{A}^1 & 0 & 0 & \mathcal{B}_{\textit{lm}}^{1} & \mathcal{B}_{\textit{lm}}^{2} & 0\\[3pt] 0 & \mathcal{B}_{il}^{3} & \mathcal{B}_{il}^{4} & 0 & 0 & \mathcal{C}_{ilmn}^{1}\\[3pt] 0 & \mathcal{B}_{il}^{4} & \mathcal{B}_{il}^{5} & 0 & 0 & \mathcal{C}_{ilmn}^{2}\\[3pt] -\mathcal{B}_{\textit{ij}}^{1} & 0 & 0 & \mathcal{C}_{\textit{ijlm}}^{3} & -\mathcal{C}_{\textit{lmij}}^{4} & \mathcal{D}_{ijlmn}^{1}\\[3pt] \mathcal{B}_{\textit{ij}}^{2} & 0 & 0 & \mathcal{C}_{\textit{ijlm}}^{4} & \mathcal{C}_{\textit{ijlm}}^{5} & \mathcal{D}_{ijlmn}^{2}\\[3pt] 0 & \mathcal{C}_{lijk}^{1} & \mathcal{C}_{\textit{ijkl}}^{2} & -\mathcal{D}_{\textit{ijklm}}^{1} & \mathcal{D}_{\textit{ijklm}}^{2} & \mathcal{E}_{ijklmn} \end{array}\right ]\left [\begin{array}{c} 1/\theta -1/T\\[3pt] (1/\theta )_{,l}\\[3pt] (1/T)_{,l}\\[3pt] v_{(l,m)}\\[3pt] - c A_{\textit{lm}}\\[3pt] - c A_{lm,n} \end{array}\right ] , \end{equation}
since the velocity gradient,
$v_{i,j}$
, and the viscous stress,
$t_{\textit{ij}}^{v}$
, are the only force and flux, respectively, that reverse under time reversal.
For isotropic gases the phenomenological coefficients are constructed from the Kronecker delta,
${\delta }_{\textit{ij}}$
, and the Levi-Civita symbol,
$\epsilon _{\textit{ijk}}$
, and can only depend on the objective scalar quantities
$\{\rho ,\theta ,T\}$
, which significantly simplifies the mathematical structure. However, the presence of the orientation tensor,
$A_{\textit{ij}}$
, introduces anisotropy that cannot be neglected. This anisotropy yields directional dependency in all fluxes, that is, thermal, granular, viscosity and orientation. It follows that the phenomenological coefficients are constructed from
$\{{\delta }_{\textit{ij}},\epsilon _{\textit{ijk}},A_{\textit{ij}}\}$
and, using the Cayley–Hamilton theorem, are functions of objective scalar state variables
$\{\rho ,\theta ,T,i_1,i_2\}$
, where
$i_1$
and
$i_2$
are the two invariants of the traceless orientational tensor,
$A_{\textit{ij}}$
. The explicit representations of the phenomenological coefficients (3.8) up to second order in
$A_{\textit{ij}}$
are given in Appendix A, which is not the most general representation but sufficient for our purpose.
4. Initial and boundary conditions
Complete solutions of the transport equations require an appropriate set of initial and boundary conditions. Initial conditions define the initial state of the system, that is, the initial values of the state variable fields throughout the domain. Boundary conditions relate the fluxes at the system’s boundaries to the properties of the boundaries and the state of the system in direct contact with them.
To construct the boundary conditions, we start with the formulation of the balance laws at a singular surface
$\mathfrak{S}$
with unit normal
$\bar {n}_i$
and velocity
$v^{\mathfrak{s}}_i$
, which separates bulk regions of different properties denoted as
$+$
and
$-$
. Specifically, we denote the granular system by
$+$
and the boundary by
$-$
. The general form of a balance law across
$\mathfrak{S}$
is expressed as a jump condition for a quantity
$\varGamma$
given by (Müller Reference Müller1985)
\begin{equation} [\![ \rho \varGamma \tilde {v}_l+\gamma _l ]\!] \bar {n}_l = \digamma {^{\mathfrak{s}}} \,, \end{equation}
where
$[\![ f ]\!]=f^{+} - f^{-}$
denotes the jump of
$f$
across
$\mathfrak{S}$
, the superscript refers to the limits at the singular surface in the direction of
$\bar {n}_i$
, the relative velocity
$\tilde {v}_i = v_i - v^{\mathfrak{s}}_i$
is the velocity of the medium with respect to
$\mathfrak{S}$
which moves with
$v^{\mathfrak{s}}_i$
,
$\gamma _i$
is the flux of the property
$\varGamma$
, and
$\digamma ^{\mathfrak{s}}$
is the production at
$\mathfrak{S}$
.
By specifying the general jump (4.1) for the variables in the balance laws, we find the jump conditions for the conservation of mass (2.3) as
\begin{equation} [\![ \rho \tilde {v}_l ]\!] \bar {n}_l =0 , \end{equation}
the conservation of linear momentum (2.4) as
\begin{equation} [\![ \rho v_i \tilde {v}_l + t_{il} ]\!] \bar {n}_l =0 , \end{equation}
the balance of thermal energy (2.7) as
\begin{equation} [\![ \rho u \tilde {v}_l + q_l ]\!] \bar {n}_l = r^{\mathfrak{s}} , \end{equation}
the balance of granular and kinetic energy (2.8) as
\begin{equation} [\![ \rho ( U + v_m v_m/2)\tilde {v}_l + Q_l + t_{\textit{ml}} v_m ]\!] \bar {n}_l = - r^{\mathfrak{s}} , \end{equation}
the balance of orientation (2.11) as
\begin{equation} [\![ A_{\textit{ij}}\tilde {v}_l + G_{\textit{ijl}} ]\!] \; \bar {n}_l =0 \end{equation}
and the balance of entropy (2.17) as
\begin{equation} [\![ \rho s \tilde {v}_l + q_l/\theta + Q_l / T -c A_{\textit{mn}} G_{mnl} ]\!] \bar {n}_l = \varSigma ^{\mathfrak{s}} \geqslant 0 , \end{equation}
where the properties of the boundary are indicated with overbar. The jump conditions for energies (4.4) and (4.5) include the exchange granular energy to thermal energy,
$r^{\mathfrak{s}}$
, across
$\mathfrak{S}$
. The entropy jump condition (4.7) includes the non-negative production term
$\varSigma ^{\mathfrak{s}}$
.
The boundary is modelled as a singular surface, hence, the jump conditions (4.2), (4.3), (4.4), (4.5) and (4.6) subject to the entropy condition (4.7), can be applied to formulate the boundary conditions. For an impermeable wall, where
$\tilde {v}_l\bar{n}_l = 0$
, the mass jump (4.2) is automatically satisfied. The remaining jump conditions reduce to
\begin{align} t_{il} \bar {n}_l &= \bar {t}_{il} \bar {n}_l , \end{align}
\begin{align} q_{l} \bar {n}_{l} &= \bar {q}_{l} \bar {n}_{l} + r^{\mathfrak{s}} , \end{align}
\begin{align} Q_{l} \bar {n}_{l} + t_{\textit{ml}}v_{m}\bar {n}_{l} &= \bar {Q}_{l} \bar {n}_{l} + \bar {t}_{\textit{ml}}\bar {v}_{m}\bar {n}_{l} - r^{\mathfrak{s}} , \end{align}
\begin{align} G_{\textit{ijl}} \bar {n}_{l} &= \bar {G}_{\textit{ijl}} \bar {n}_{l} . \end{align}
Substitution of these relations to eliminate the boundary fluxes
$\{\bar {t}_{\textit{ij}},\bar {q}_i,\bar {Q}_i,\bar {G}_{\textit{ijk}}\}$
in the entropy jump (4.7) yields
\begin{align} \varSigma ^{\mathfrak{s}} &= q_l \bar {n}_l ( 1 / \theta - 1 / \bar {\theta } ) + r^{\mathfrak{s}} (1/\bar {\theta }-1/\bar {T}) + Q_l \bar {n}_l ( 1 / T - 1 / \bar {T} ) \nonumber \\&\!\quad + t^\tau _{\textit{ml}}\bar {n}_l/\bar {T} \tilde {v}_m^\tau + G_{mnl} \bar {n}_l c ( \bar {A}_{\textit{mn}} -A_{\textit{mn}} ) \geqslant 0 , \end{align}
where
$( 1 / \theta - 1 / \bar {\theta } )$
is the temperature jump,
$( 1 / T - 1 / \bar {T} )$
is the granular temperature jump and
$( A_{\textit{ij}} - \bar {A}_{\textit{ij}} )$
is the orientation jump between the grains and the wall. The tangential velocity jump, also called the tangential slip velocity, is
$\tilde {v}_i^\tau = v_i - \bar {v}_i$
where
$\tilde {v}_l^\tau \bar {n}_l=0$
. Interestingly, the term
$(1/\bar {\theta }-1/\bar {T})$
does not represent a jump, but rather indicates the difference between the thermal and granular temperatures at the wall, governing the exchange of granular energy into thermal energy, denoted by
$r^{\mathfrak{s}}$
. The fluxes into the system include heat flux,
$q_l \bar {n}_l$
, granular heat flux,
$Q_l \bar {n}_l$
, linear momentum flux,
$t^{\tau }_{il}\bar {n}_{l}$
, and orientation flux,
$G_{\textit{ijl}} \bar {n}_l$
.
The wall temperature,
$\bar {\theta }$
, is a well-defined property, but the granular temperature,
$\bar {T}$
, and the orientation,
$\bar {A}_{\textit{ij}}$
, of the wall need to be defined. By considering the wall as being constructed of the same type of particles as those in the system, the effective granular temperature and orientation of the wall can be defined as those of the particles composing it.
By the second law of thermodynamics, the non-negative entropy production across the singular surface, given in (4.12), is the sum of products of boundary forces and fluxes
\begin{align} \varSigma ^{\mathfrak{s}}=\sum _{\alpha } \mathcal{F}^{\mathfrak{s}}_\alpha \mathcal{J}^{\mathfrak{s}}_\alpha \geqslant 0 \,. \end{align}
From (4.12) , we identify the forces and fluxes as
\begin{align} \mathcal{F}^{\mathfrak{s}}_\alpha &= \begin{bmatrix} \; 1/\theta - 1/\bar {\theta } , 1/\bar {\theta }-1/\bar {T} , 1/T - 1/\bar {T} , \tilde {v}^\tau _i ,\; c(\bar {A}_{\textit{ij}} - A_{\textit{ij}}) \; \end{bmatrix} \!, \end{align}
\begin{align} \mathcal{J}_\alpha ^{\mathfrak{s}} &= \begin{bmatrix} \; q_l \bar {n}_l , r^{\mathfrak{s}} , Q_l \bar {n}_l , t^\tau _{il} \bar {n}_l / \bar {T} , G_{\textit{ijl}} \bar {n}_l \; \end{bmatrix} \!, \end{align}
respectively, where the forces vanish at equilibrium.
Since the boundary fluxes
$\mathcal{J}_\alpha ^{\mathfrak{s}}$
associated with irreversible processes (De Groot & Mazur Reference De Groot and Mazur1962) must be proportional to the boundary forces,
$\mathcal{F}_\alpha ^{\mathfrak{s}}$
, they are taken as linear functions of the boundary forces
\begin{align} \mathcal{J}_{\alpha }^{\mathfrak{s}}=\sum _{\beta } L_{\alpha \beta }^{\mathfrak{s}} \mathcal{F}_\beta ^{\mathfrak{s}} \,, \end{align}
where
$L^{\mathfrak{s}}_{\alpha \beta }$
is the matrix of phenomenological coefficients that characterizes the constitutive properties of the boundary interactions. To ensure non-negative entropy production, consistent with (3.6), the matrix
$L^{\mathfrak{s}}_{\alpha \beta }$
must be positive semidefinite. Invoking the Onsager–Casimir reciprocity principle (Onsager Reference Onsager1931a
,
Reference Onsagerb
; Casimir Reference Casimir1945), the force-flux relations are
\begin{equation} \begin{bmatrix} q_l \bar {n}_l \\[3pt] r^{\mathfrak{s}} \\[3pt] Q_l \bar {n}_l \\[3pt] t^\tau _{il} \bar {n}_l/\bar {T} \\[3pt] G_{\textit{ijl}} \bar {n}_l \end{bmatrix} = \begin{bmatrix} \mathcal{A}^2 & \mathcal{A}^3 & \mathcal{A}^4 & \mathcal{H}^1_{l} & \mathcal{B}^6_{\textit{lm}} \\[3pt] \mathcal{A}^3 & \mathcal{A}^5 & \mathcal{A}^6 & \mathcal{H}^2_{l} & \mathcal{B}^7_{\textit{lm}} \\[3pt] \mathcal{A}^4 & \mathcal{A}^6 & \mathcal{A}^7 & \mathcal{H}^3_{l} & \mathcal{B}^8_{\textit{lm}} \\[3pt] -\mathcal{H}^1_{i} & -\mathcal{H}^2_{i} & -\mathcal{H}^3_{i} & \mathcal{B}^9_{il} & \mathcal{G}_{ilm} \\[3pt] \mathcal{B}^6_{\textit{ij}} & \mathcal{B}^7_{\textit{ij}} & \mathcal{B}^8_{\textit{ij}} & -\mathcal{G}_{\textit{ijl}} & \mathcal{C}^6_{\textit{ijlm}} \end{bmatrix} \begin{bmatrix} 1/\theta - 1/\bar {\theta } \\[3pt] 1/\bar {\theta } - 1/\bar {T} \\[3pt] 1/T - 1/\bar {T} \\[3pt] \tilde {v}_l^\tau \\[3pt] \bar {A}_{\textit{lm}} - A_{\textit{lm}} \end{bmatrix} . \end{equation}
The phenomenological coefficients in
$L^{\mathfrak{s}}_{\alpha \beta }$
depend on the objective state variables of the system
$\{\rho ,\theta ,T,A_{\textit{ij}}\}$
and on the properties of the boundary
$\{\bar {\theta },\bar {T},\bar {n}_i,\bar {A}_{\textit{ij}}\}$
. As such, the phenomenological coefficients are constructed using (Wang Reference Wang1969) the following additional tensors
$\{\bar {n}_i,\bar {A}_{\textit{ij}}\}$
and the additional objective scalar invariants
\begin{align} &i_3 = A_{\textit{lm}}\bar {n}_l \bar {n}_m , i_4 = A^2_{\textit{lm}} \bar {n}_l \bar {n}_m , i_5 = A_{\textit{lm}}\bar {A}_{\textit{lm}} , i_6 = A_{\textit{lm}}\bar {A}^2_{\textit{lm}} , i_7 = A^2_{\textit{lm}}\bar {A}_{\textit{lm}} , \nonumber \\& i_8 = A^2_{\textit{lm}}\bar {A}^2_{\textit{lm}} , i_9 = \bar {A}_{\textit{lm}} \bar {n}_l \bar {n}_m , i_{10} = \bar {A}^2_{ll} , i_{11} = \bar {A}^2_{\textit{lm}} \bar {n}_l \bar {n}_m , i_{12} = \bar {A}^3_{ll} , \end{align}
where
$\{i_3,i_4,i_5,i_6,i_7,i_8\}$
are additional invariants of the orientation field at the boundary and
$\{i_9,i_{10},i_{11},i_{12}\}$
are invariant properties of the boundary. Hence, the objective scalar state variables and boundary properties are
\begin{equation} \bigl\{\rho ,\theta ,T,\bar {\theta },\bar {T},i_{1-12}\bigr\} . \end{equation}
It is worth noting that the tacit consideration in taking the fluxes to be linear functions of the forces in (4.17) is that the non-equilibrium state is a small deviation from the equilibrium state, where,
$\theta \approx \bar {\theta }$
,
$T \approx \bar {T}$
and
$A_{\textit{ij}} \approx \bar {A}_{\textit{ij}}$
, hence, the boundary properties
$\{\bar {\theta },\bar {T},\bar {A}_{\textit{ij}}\}$
are typically not considered explicitly in the phenomenological coefficients taken at equilibrium state. However, we do not make this assumption here, allowing the formulation in (4.17) to be, potentially, applicable farther from equilibrium.
5. Applications to oriented granular gases
The general constitutive laws (3.8) are limited by taking the fluxes (3.3) to be linear functions of the forces (3.2) but otherwise general functions of the objective state variables (3.5). The phenomenological coefficients for the bulk and the boundaries as presented in Appendices A and B involve a total of
$292$
transport coefficients to be determined. Clearly, this level of generality, while theoretically comprehensive, is impractical for direct application and is best regarded as a reference framework for future generalization of the constitutive laws, if needed.
Next, we demonstrate that available results in the literature of granular materials (Jenkins Reference Jenkins2007; Vescovi et al. Reference Vescovi, Nadler and Berzi2024; Berzi et al. Reference Berzi, Vescovi and Nadler2025) are simplified cases of the general framework developed here. To obtain a simpler and physically interpretable model, we introduce significant simplifications by systematically neglecting what we consider secondary effects. In particular, we eliminate most of the off-diagonal coupling terms in the force-flux relations, retaining only the dominant contributions
\begin{equation} \left [\begin{array}{c} r\\ q_{i}\\[3pt] Q_{i}\\[3pt] t_{\textit{ij}}^{v}/T\\[3pt] P_{\textit{ij}}\\[3pt] G_{\textit{ijk}} \end{array}\right ]=\left [\begin{array}{cccccc} \mathcal{A}^1 & 0 & 0 & 0 & 0 & 0\\[3pt] 0 & \mathcal{B}_{il}^{3} & 0 & 0 & 0 & 0 \\[3pt] 0 & 0 & \mathcal{B}_{il}^{5} & 0 & 0 & 0 \\[3pt] 0 & 0 & 0 & \mathcal{C}_{\textit{ijlm}}^{3} & -\mathcal{C}_{\textit{lmij}}^{4} & 0 \\[3pt] 0 & 0 & 0 & \mathcal{C}_{\textit{ijlm}}^{4} & \mathcal{C}_{\textit{ijlm}}^{5} & 0\\[3pt] 0 & 0 & 0 & 0 & 0 & \mathcal{E}_{ijklmn} \end{array}\right ]\left [\begin{array}{c} 1/\theta -1/T\\[3pt] (1/\theta )_{,l}\\[3pt] (1/T)_{,l}\\[3pt] v_{(l,m)}\\[3pt] -cA_{\textit{lm}}\\[3pt] -cA_{lm,n} \end{array}\right ] , \end{equation}
and similarly for the boundaries
\begin{equation} \begin{bmatrix} q_l \bar {n}_l \\[3pt] r^{\mathfrak{s}} \\[3pt] Q_l \bar {n}_l \\[3pt] t^\tau _{il} \bar {n}_l/\bar {T} \\[3pt] G_{\textit{ijl}} \bar {n}_l \end{bmatrix} = \begin{bmatrix} \mathcal{A}^2 & 0 & 0 & 0 & 0 \\[3pt] 0 & \mathcal{A}^5 & 0 & 0 & 0 \\[3pt] 0 & 0 & \mathcal{A}^7 & 0 & 0 \\[3pt] 0 & 0 & 0 & \mathcal{B}^9_{il} & 0 \\[3pt] 0 & 0 & 0 & 0 & \mathcal{C}^6_{\textit{ijlm}} \end{bmatrix} \begin{bmatrix} 1/\theta - 1/\bar {\theta } \\[3pt] 1/\bar {\theta } - 1/\bar {T} \\[3pt] 1/T - 1/\bar {T} \\[3pt] \tilde {v}^\tau _l \\[3pt] \bar {A}_{\textit{lm}} - A_{\textit{lm}} \end{bmatrix} . \end{equation}
This simplification reduces the number of transport coefficients to
$153$
, hence more simplifications are required. Simplification is obtained by comparison with results available in the literature for spherical and non-spherical particles to identify some of the transport coefficients. In the absence of such results we propose simple generalization by including dependency of the coefficients on the orientation.
The second law of thermodynamic requires that the tensors
$\{\mathcal{A}^1,\mathcal{B}_{\textit{lm}}^{3},\mathcal{B}_{\textit{lm}}^{5},\mathcal{C}_{\textit{ijlm}}^{3}, \mathcal{C}_{\textit{ijlm}}^{5},\mathcal{E}_{ijklmn}\}$
and
$\{\mathcal{A}^2,\mathcal{A}^5,\mathcal{B}_{\textit{lm}}^{9},\mathcal{C}_{\textit{ijlm}}^{6}\}$
be positive semidefinite imposing constraints on the transport coefficients, however, no such restriction applies to the tensor
$\mathcal{C}_{\textit{ijlm}}^{4}$
. The phenomenological coefficients must depend on the objective scalar state variables
$\{\rho ,\theta ,T,i_1,i_2\}$
, and at the boundaries they must also depend on the coupling invariants
$\{i_{3-6}\}$
and on the boundary properties
$\{\bar {\theta },\bar {T},i_{7-12}\}$
. In addition, the bounds of the orientation tensor provide further restrictions, (2.15) and (2.16), on the admissible range of the transport coefficients.
5.1. Determination of the transport coefficients of oriented granular gases
In this section the transport coefficients are determined by comparison with results available in the literature for spherical and non-spherical particles (Soto Reference Soto2004; Jenkins Reference Jenkins2007; Vescovi et al. Reference Vescovi, Nadler and Berzi2024), as well as by proposing simple extensions when such results are not available. For granular gases, the relevant particle properties include geometry, mass and the response to collisions. In the granular kinetic theory of spherical particles, the key properties are
$\{\rho _p,d_p,e_p\}$
, where
$\rho _p$
is the solid particle density,
$d_p$
is the particle diameter and
$e_p$
is the coefficient of restitution. Collisions are elastic for
$e_p = 1$
, inelastic for
$e_p \in [0,1)$
and perfectly plastic for
$e_p = 0$
. For oriented particles, the deviation from spherical shape can be measured by the aspect ratio
$r_p=(l-d)/(l+d)$
, where
$l$
is the length along the particle axis and
$d$
is the diameter. Hence,
$r_p\in [-1,1]$
, where
$r_p\gt 0$
are prolate particles and
$r_p\lt 0$
are oblate particles. The equivalent particle diameter
$d_p$
can be defined as the diameter of a spherical particles of similar volume. It should be noted that, in granular gases, elastic properties do not play any role as particle interactions are dominated by collisions rather than long-term contacts. Similarly, the boundary property
$e_w$
is the coefficient of restitution between the particles and the boundary.
While an extensive body of work exists for spherical particles, studies addressing oriented particles remain limited. By taking spherical particles as a special case, where
$A_{\textit{ij}}=0$
and
$r_p=0$
, our general formulation, incorporating the orientation, should reduce to the established results for spherical particles while allowing generalization to oriented particles. Motivated by findings from the granular kinetic theory for spherical particles, where
$A_{\textit{ij}}=0$
, such as those presented in Jenkins (Reference Jenkins2007) and by recent extension to oriented particles presented in Vescovi et al. (Reference Vescovi, Nadler and Berzi2024), we propose forms for some of the phenomenological transport coefficients.
For granular materials, it is convenient to define the dimensionless solid volume fraction by
\begin{equation} \varrho = \rho / \rho _p , \end{equation}
as it plays an important role.
5.2. Transport coefficients of granular gases of oriented particles
The transport coefficients from granular gases of oriented particles (5.1) are proposed by generalization of available results for the granular kinetic theory of spherical particles and other available results for oriented particles. Finally, a summary of the system of transport equations is provided in § 6.
5.2.1. Thermal and granular energies
The granular kinetic theory of spherical particles provides explicit expressions for the transport coefficients associated with granular energy. However, the thermal transport coefficients are not addressed explicitly since the thermal energy is considered to be much smaller than the granular energy.
However, we retain an explicit consideration of the thermal energy in order to provide a more general and self-consistent thermodynamic framework. This allows the theory to be applied to situations where thermal effects, such as heat transfer, are relevant, and it also permits the transport coefficients to depend on the thermal temperature, for example through a temperature-dependent coefficient of restitution. In the limit where thermal effects are negligible, the formulation naturally reduces to the commonly used granular models.
The exchange of granular energy into thermal energy due to inelastic collisions (5.1) is
\begin{equation} r = \mathcal{A}^1 (1/\theta -1/T) . \end{equation}
When taking
$T \gg \theta$
, which is relevant to most granular gases, and adopting a result from the granular kinetic theory of spherical particles (Jenkins Reference Jenkins2007),
$A_{\textit{ij}}=0$
, it simplifies to
\begin{equation} r = \mathcal{A}^1 (1/\theta ) = \frac {24 (1-e_p)}{\sqrt {\pi }d_p} \rho \chi \sqrt {T^3} , \end{equation}
where
$\chi (\varrho )$
is the radial distribution function given in Jenkins (Reference Jenkins2007). We propose a simple generalization of (5.5) to oriented particles by introducing a scalar correction factor
$\alpha _1\geqslant 0$
which depends on the orientation and aspect ratio as
\begin{equation} \mathcal{A}^1 = \alpha _1(i_1,i_2,r_p)\frac {24 (1-e_p)}{\sqrt {\pi }d_p}\rho \chi \sqrt {T^3} \theta , \alpha _1\geqslant 0 , \alpha _1(0) = 1 , \end{equation}
which satisfies the non-negative requirement. The precise form of
$\alpha _1(i_1,i_2,r_p)$
remains an open question for future investigation.
The granular heat flux (5.1) is given by
\begin{equation} Q_i = \mathcal{B}^5_{il} (1/T)_{,l} , \end{equation}
where the granular heat conductivity tensor
$\mathcal{B}^5_{\textit{ij}}$
is taken to be linear in the orientation (A3)
\begin{equation} \mathcal{B}^{5}_{\textit{ij}} = b_5 {\delta }_{\textit{ij}} + \check {b}_5 A_{\textit{ij}} , b_5 \gt 0 , -3/2 \lt \check {b}_5/b_5 \lt 3 . \end{equation}
We adopt a result from granular the kinetic theory of spherical particles (Jenkins Reference Jenkins2007),
$A_{\textit{ij}}=0$
, as
\begin{equation} Q_i = \mathcal{B}^5_{il} (1/T)_{,l} = \frac {4 d_p}{\sqrt {\pi }}\rho M \chi \sqrt {T} T_{,i} \end{equation}
where the functional form of
$M$
is given in Jenkins (Reference Jenkins2007). We propose a simple generalization of (5.9) to oriented particles by introducing a scalar correction factor
$\beta _5\geqslant 0$
which depends on the orientation and aspect ratio. It follows that
\begin{equation} b_5 = \beta _5(i_1,i_2,r_p)\frac {4 d_p}{\sqrt {\pi }}\rho M \chi \sqrt {T^3} , \beta _5 \geqslant 0 ,\beta _5(0)=1 , \end{equation}
where the explicit forms of
$\beta _1(i_1,i_2,r_p)$
and
$\check {b}_5$
remain open questions for future investigation.
Next, the thermal heat flux is given by
\begin{equation} q_i = \mathcal{B}^3_{il} (1/\theta )_{,l} . \end{equation}
Thermal heat conduction typically plays a minor role in granular gases, as
$T\gg \theta$
, and it is not addressed in the granular kinetic theory, hence, this phenomenological coefficient is not available from granular gases. However, taking a linear dependence on the orientation tensor (A3), the thermal heat conductivity is
\begin{equation} \mathcal{B}^{3}_{\textit{ij}} = b_3 {\delta }_{\textit{ij}} + \check {b}_3 A_{\textit{ij}} , b_3 \gt 0 , -3/2 \lt \check {b}_3/b_3 \lt 3 , \end{equation}
where determining the transport coefficients
$b_3$
and
$\check {b}_3$
as functions of the objective scalar state variables (3.5) requires further investigation.
The transport coefficients
$\check {b}_3$
and
$\check {b}_5$
represent the anisotropy in the thermal and granular heat conductivities, respectively, induced by the orientation tensor. While no qualitative or quantitative data is currently available, it is intuitive that particle alignment can induce anisotropic conductivities, i.e. positive values of these coefficients increase conductivity along the alignment direction and reduce it perpendicular to the alignment. Conversely, negative values of these coefficients reduce conductivity in the alignment direction and increased conductivity perpendicular to it.
5.2.2. Stress and orientation production
The phenomenological coefficients in (5.1), which govern the coupled viscous stress and orientation production, are fourth-order tensors. The viscous stress is given by
\begin{equation} t^v_{\textit{ij}}/T = \mathcal{C}^3_{\textit{ijlm}} v_{(lm)} + c \mathcal{C}^4_{\textit{lmij}} A_{\textit{lm}} , \end{equation}
and the orientation production is
\begin{equation} P_{\textit{ij}} = \mathcal{C}^4_{\textit{ijlm}} v_{(lm)} - c \mathcal{C}^5_{\textit{ijlm}} A_{\textit{lm}} . \end{equation}
Following Vescovi et al. (Reference Vescovi, Nadler and Berzi2024), the forms of
$\mathcal{C}^3_{\textit{ijlm}}$
,
$\mathcal{C}^4_{\textit{ijlm}}$
and
$\mathcal{C}^5_{\textit{ijlm}}$
are determined for cylindrical particles and no further generalizations are proposed here.
The phenomenological coefficient
$\mathcal{C}^3_{\textit{ijkl}}$
given in (A5) is expressed using five terms
\begin{equation} \mathcal{C}^3_{\textit{ijlm}} = c_{3_1} {\delta }_{\textit{ij}}{\delta }_{\textit{lm}} + c_{3_2} {\delta }_{(i(l}{\delta }_{m)j)} + \check {c}_{3_1} {\delta }_{\textit{ij}}A_{\textit{lm}} + \check {c}_{3_2} A_{\textit{ij}}{\delta }_{\textit{lm}} + \check {c}_{3_3} {\delta }_{(i(l}A_{m)j)} , \end{equation}
where the coefficients are related to three transport parameters
$\{\lambda ,\mu ,\eta \}$
as
\begin{equation} c_{3_1} = \lambda /T , c_{3_2} = 2\mu /T , \check {c}_{3_1} = -2\mu \eta /(3T) , \check {c}_{3_2} = -2 \mu \eta /T , \check {c}_{3_3} = 3 \mu \eta /T . \end{equation}
Here,
$\lambda$
and
$\mu$
are the bulk and shear viscosities, respectively, and
$\eta$
represents the anisotropic effect of orientation on the shear viscosity. The positive semidefinite requirement for
$\mathcal{C}^3_{\textit{ijkl}}$
is satisfied for
\begin{equation} \lambda \geqslant 0 , \mu \geqslant 0 , \eta \in [0,1] , \end{equation}
hence,
$\eta$
represents the reduction in shear viscosity along the alignment direction. As shown in Vescovi et al. (Reference Vescovi, Nadler and Berzi2024), for frictionless cylinders,
$\eta =1$
, and it decreases with increasing surface roughness of the particles.
The bulk and shear viscosities are adopted from the granular kinetic theory of spherical particles (Jenkins Reference Jenkins2007)
\begin{equation} \lambda = \frac {4 d_p (1+e_p)}{3\sqrt {\pi }}\rho \chi \sqrt {T} , \mu = \frac {8 d_p}{5\sqrt {\pi }}\rho \chi \sqrt {T} . \end{equation}
The pressure
$p$
is found to be independent of the orientation and aspect ratio and in good agreement with the granular kinetic theory of spherical particles (Berzi et al. Reference Berzi, Thai-Quang, Guo and Curtis2016), and takes (Jenkins Reference Jenkins2007) the form
\begin{equation} p = \rho ( 1 + 4 \chi ) T . \end{equation}
It is somewhat expected that pressure is independent of orientation due to its isotropic nature.
It is important to note that the term
$c\mathcal{C}^4_{\textit{lmij}} A_{\textit{lm}}$
in (5.13) is neglected in Vescovi et al. (Reference Vescovi, Nadler and Berzi2024) when considering steady-state shear flows. While further investigation is needed, it is expected that for sufficiently fast shear flow, where
$\mathcal{C}^3_{\textit{ijlm}} v_{(l,m)} \gg c\mathcal{C}^4_{\textit{lmij}} A_{\textit{lm}}$
, the stress response is dominated by the term
$\mathcal{C}^3_{\textit{ijlm}} v_{(l,m)}$
. However, the role of the term
$c\mathcal{C}^4_{\textit{lmij}} A_{\textit{lm}}$
becomes important as
$v_{(l,m)}\rightarrow 0$
, such as granular agitation in confined volumes.
We proceed by identifying the phenomenological coefficients
$\mathcal{C}^4_{\textit{ijkl}}$
and
$\mathcal{C}^5_{\textit{ijkl}}$
governing the orientation production. In Vescovi et al. (Reference Vescovi, Nadler and Berzi2024),
$\mathcal{C}^4_{\textit{ijkl}}$
is expressed as a combination of four terms (A6) involving a single transport coefficient
\begin{equation} \mathcal{C}^4_{\textit{ijlm}} = c_{4} {\delta }_{\langle i(l}{\delta }_{m)j\rangle } + \check {c}_{4_1} A_{\textit{ij}}{\delta }_{\textit{lm}} + \check {c}_{4_2} A_{\langle i(l}{\delta }_{m)j\rangle } + \hat {c}_{4_1} A_{\textit{ij}}A_{\textit{lm}} , \end{equation}
with the coefficients taken to be
\begin{equation} c_{4} = 2/3 \rho \phi , \check {c}_{4_1} = -2/3 \rho \phi , \check {c}_{4_2} = 2 \rho \phi , \hat {c}_{4_1} = -2 \rho \phi . \end{equation}
It should be noted that the non-negative entropy production in (5.1) imposes no restrictions on
$\phi (r_p)$
, however, as argued in Nadler et al. (Reference Nadler, Guillard and Einav2018) and confirmed via discrete simulations in Vescovi et al. (Reference Vescovi, Nadler and Berzi2024),
$\phi (r_p)\in [-1,1]$
and
$\phi (-r_p)\approx -\phi (r_p)$
.
The phenomenological coefficient
$\mathcal{C}^5_{\textit{ijlm}}$
is expressed (A7) using a single term and a single transport coefficient
\begin{equation} c\mathcal{C}^5_{\textit{ijlm}} = c_{5} {\delta }_{\langle i\langle l}{\delta }_{m\rangle j\rangle } , \end{equation}
where the positive semidefinite requirement of
$c\mathcal{C}^5_{\textit{ijlm}}$
is satisfied for
$c_{5} \geqslant 0$
. Furthermore, Vescovi et al. (Reference Vescovi, Nadler and Berzi2024) proposes the following form:
\begin{equation} c_{5}=\psi (\varrho ,r_p) d_p^{-1}\rho \sqrt {T} , \; \psi (\varrho ,r_p) \geqslant 0 , \end{equation}
where
$\psi (\varrho ,r_p)$
is of order one, and
$\psi (\varrho ,-r_p)\approx \psi (\varrho ,r_p)$
.
It is shown in Amereh & Nadler (Reference Amereh and Nadler2023) that for (5.20), (5.21) and (5.22), the restriction (2.15) is satisfied, where the coefficient
$\hat {c}_{4_1}$
corresponding to the quadratic term in (5.20) is essential to enforce orientation bounds in the case of highly aligned particles.
5.2.3. Orientation flux
By (5.1) the orientation flux is
\begin{equation} G_{\textit{ijk}} = -c\mathcal{E}_{ijklmn} A_{lm,n} . \end{equation}
In the absence of meaningful data to guide the form of these coefficients, we adopt the simplest form, where the orientation flux
$G_{\textit{ijk}}$
is aligned with the orientation gradient and depends on a single transport coefficient. By inspecting (A10), we retain only one term,
$e_3=\rho k$
, leading to
\begin{equation} c\mathcal{E}_{ijklmn} = \rho k {\delta }_{\langle i\langle l}{\delta }_{m\rangle j\rangle }{\delta }_{kn} , k \gt 0 . \end{equation}
More specifically, since transport coefficients in the granular kinetic theory (Jenkins Reference Jenkins2007) are proportional to the collision frequency, we propose
\begin{equation} k = \kappa (\varrho , i_1,i_2,r_p) d_p \sqrt {T} , \kappa \geqslant 0 , \end{equation}
where
$k$
is proportional to
$\sqrt {T}$
, similar to (5.10) and (5.18). However, no explicit form for the dimensionless coefficient
$\kappa (\varrho , i_1,i_2,r_p)$
is currently available. The positive semidefinite requirement of
$c\mathcal{E}_{ijklmn}$
is satisfied, and it is shown in Amereh & Nadler (Reference Amereh and Nadler2023) that also (2.16) is fulfilled.
5.3. Boundary conditions of oriented granular gases
The boundary conditions for granular gases of oriented particles (5.2) are proposed by generalizing results from the granular kinetic theory of spherical particles, when available, to describe the exchange of granular energy into thermal energy due to inelastic interactions with the boundary, granular heat flux, thermal heat flux, momentum flux and orientation flux from the boundary.
5.3.1. Boundary energy exchange and heat fluxes
The exchange of granular energy into thermal energy due to interactions with the boundary (5.2) is given by
\begin{equation} r^{\mathfrak{s}} = \mathcal{A}^{5}( 1/\bar {\theta } - 1/\bar {T} ) , \end{equation}
where the non-negative requirement is satisfied for
$\mathcal{A}^{5}\geqslant 0$
. Equation (5.27) indicates that due to the interaction of the particles with the wall, when
$\bar {T}\gt \bar {\theta }$
granular energy transfers to thermal energy, when
$\bar {\theta }\gt \bar {T}$
thermal energy transfers to granular energy, and when
$\bar {\theta }=\bar {T}$
the transfer between granular energy to thermal energy vanishes. In the relevant limit where
$\bar {T} \gg \bar {\theta }$
, a result from (Jenkins Reference Jenkins2007) the granular kinetic theory of spherical particles,
$A_{\textit{ij}}=0$
, is
\begin{equation} r^{\mathfrak{s}} = \mathcal{A}^{5}( 1/\bar {\theta } ) = \sqrt {2/\pi }(1-e_w)\rho (1+4\chi ) \sqrt {T^3} , \end{equation}
where
$e_w$
is the coefficient of restitution between a particle and the wall. We generalized it to oriented particles by multiplication with a correction factor that accounts for the orientation and aspect ratio yielding
\begin{equation} \mathcal{A}^{5} = \alpha _5(i_{1-8},r_p) \sqrt {2/\pi }(1-e_w)\rho (1+4\chi ) \sqrt {T^3}\bar {\theta } , \alpha _5 \geqslant 0 ,\; \alpha _5(0) = 1 . \end{equation}
The form of
$\alpha _5(i_{1-8},r_p)$
is an open question that requires further investigation.
The granular heat flux at the boundary is
\begin{equation} Q_l \bar {n}_l = \mathcal{A}^{7} (1/T - 1/\bar {T}) , \end{equation}
which is proportional to the granular temperature jump. In Soto (Reference Soto2004), the granular kinetic theory is used to evaluate this flux for a vibrating wall in contact with spherical particles, in term of the amplitude
$\mathbb{A}$
and the frequency
$w$
of vibrations, in the limit of small amplitude
$\mathbb{A}\rightarrow 0$
and high frequency
$w\rightarrow \infty$
such that the product
$\mathbb{A}w$
remains finite. The granular temperature of the wall is
$\bar {T}=(\mathbb{A}w)^2$
and the inward granular heat flux from the wall becomes
\begin{equation} Q_l \bar {n}_l = \frac {4}{5} \rho \bar {T} \sqrt {T} , \end{equation}
which does not explicitly account for the granular temperature jump. To remain consistent with the second law of thermodynamic where, depends on the granular temperature jump, allowing both inward and outward fluxes, it is considered in Soto (Reference Soto2004) that
$\bar {T} \gg T$
so the granular heat flux is directed into the granular system, hence
\begin{equation} Q_l \bar {n}_l = \mathcal{A}^{7} (1/T) = \frac {4}{5} \rho \bar {T} \sqrt {T} , \end{equation}
which is valid for spherical particles. Hence, we propose a generalization to oriented particles as
\begin{equation} \mathcal{A}^{7} = \alpha _7(i_{1-8},r_p) \frac {4}{5} \rho \bar {T} \sqrt {T^3} , \alpha _7 \geqslant 0 , \alpha _7(0) =1 , \end{equation}
where the form of
$\alpha _7(i_{1-8},r_p)$
remains an open question. However, due to the isotropic nature of this effect, a strong dependence on orientation is not expected, instead, the dependence is expected to be primarily on particle shape
$\alpha _7(r_p)$
, but this must be investigated.
The boundary thermal heat flux governed by the phenomenological coefficient
$\mathcal{A}^{2}$
is not addressed in the granular kinetic theory and generally plays a minor role in granular gases. This transport coefficient can be a function of all objective state variables and the boundaries, subject to the requirement
\begin{equation} \mathcal{A}^{2}(\rho ,\theta ,T,\bar {\theta },\bar {T},i_{1-8}) \geqslant 0 , \end{equation}
however, no explicit expression is currently available.
5.3.2. Boundary traction and slip velocity
The relationship between the tangential boundary traction and the boundary slip velocity (5.2) is
\begin{equation} t^\tau _{il} \bar {n}_l/\bar {T} = \mathcal{B}^9_{il} \tilde {v}^\tau _l . \end{equation}
Taking (B4) to be linear in the orientation yields
\begin{align} \mathcal{B}^{9}_{\textit{ij}} = b_9 {\delta }_{\textit{ij}} + \check {b}_9 A_{\textit{ij}} , b_9 \gt 0 , -3/2 \lt \check {b}_9/b_9 \lt 3 . \end{align}
It is reasonable to consider that the orientation of both the boundary and the particles contribute equally, which can be simply incorporated by
\begin{align} \mathcal{B}^{9}_{\textit{ij}} = b_9 {\delta }_{\textit{ij}} + \frac {\check {b}_9}{2} \bigl( A_{\textit{ij}} + \bar {A}_{\textit{ij}}\bigr) , b_9 \gt 0 , -3/4 \lt \check {b}_9/b_9 \lt 3/2 . \end{align}
However, since the coefficient
$\check {b}_9$
is unknown and we are considering small deviations from equilibrium where
$A_{\textit{ij}} \approx \bar {A}_{\textit{ij}}$
, it is sufficient to use (5.36). We adopt a result from (Jenkins Reference Jenkins2007) the granular kinetic theory of spherical particles,
$A_{\textit{ij}}=0$
, where
\begin{equation} t^\tau _{il} \bar {n}_l/\bar {T} = \mathcal{B}^9_{il} \tilde {v}^\tau _l = \frac {\zeta ^2}{\sqrt {2\pi }} (1+4\chi ) \rho \sqrt {T}\tilde {v}^\tau _i , \end{equation}
where, the boundary bumpiness
$\zeta$
measures the average maximum penetration of flow particles into the boundary asperities. It follows that a generalization to oriented particles is
\begin{equation} b_9 = \beta _{9}(i_{1-8},r_p)\frac {\zeta ^2}{\sqrt {2\pi }} (1+4\chi ) \rho \sqrt {T} , \beta _{9} \geqslant 0 , \beta _{9}(0) = 1 , \end{equation}
where the forms of the transport coefficients
$\beta _9$
and
$\check {b}_9$
remain open questions and require further investigation.
5.3.3. Boundary orientation flux
The orientation flux from the boundary (5.2) is given by
\begin{equation} G_{\textit{ijl}} \bar {n}_l = \mathcal{C}^6_{\textit{ijlm}} (\bar {A}_{\textit{lm}} - A_{\textit{lm}} ) , \end{equation}
suggesting that particles tend to align with the boundary. Here we consider the orientation of a wall as the orientation of the particles constructing it, where these particles are taken to be identical to the oriented particles of the system.
This provides a clear definition for the orientation of a flat wall when the system consists of oblate particles, hence
$\bar {A}_{\textit{ij}}=\bar {n}_i \bar {n}_j - {\delta }_{\textit{ij}}/3$
. However, this definition becomes ambiguous for prolate particles since the boundary is formed by in-plane particles whose specific in-plane orientations are undefined. As such,
$\bar {A}_{\textit{ij}}$
is not fully defined for prolate particles, but rather, the only know property is
$\bar {A}_{il} \bar {n}_l= -1/3 \bar {n}_i$
, which indicates that no particles at the wall are oriented normal to it.
In the absence of any definitive data, similar to (5.25), we adopt the simplest form, where the orientation flux is proportional and aligned with the orientation jump. From (B7), this yields
\begin{equation} \mathcal{C}^6_{\textit{ijlm}} = c_{6_1} {\delta }_{\langle i\langle l}{\delta }_{m\rangle j\rangle } + c_{6_2} \bar {n}_{\langle i}{\delta }_{j\rangle \langle l}\bar {n}_{m\rangle } . \end{equation}
For oblate particles in contact with a flat wall, this is
\begin{equation} c_{6_1} = \rho k_w , c_{6_2}=0 , \bar {A}_{\textit{ij}} = \bar {n}_i \bar {n}_j - \frac {1}{3}{\delta }_{\textit{ij}} , \end{equation}
while for prolate particles it is
\begin{equation} c_{6_1} = 0 , c_{6_2} = \rho k_w , \bar {A}_{il} \bar {n}_l= -1/3 \bar {n}_i . \end{equation}
Following the same argument as in (5.26), we take the transport coefficient as
\begin{equation} k_w = \kappa _w(\varrho ,i_{1-8},r_g) d_p^{-1} \sqrt {T} , \kappa _w \geqslant 0 . \end{equation}
The form of the coefficient
$\kappa _w$
is an open question requiring further investigation.
6. Summary of the transport equations of oriented granular gases
We derived a set of the governing equations for anisotropic granular gases that can be used to solve initial-boundary value problems. In granular gases, it is typically taken that
$T\gg \theta$
making the thermal energy negligible, thus, it is not explicitly considered. The conservation and balance laws for the state variables
$\{\rho ,v_i,T,A_{\textit{ij}}\}$
are
\begin{align} \dot {\rho } + \rho v_{i,i} &=0 , \end{align}
\begin{align} \rho \dot {v}_i + t_{\textit{ij},j} &= 0 , \end{align}
\begin{align} \rho \frac {5}{2}\dot {T} + Q_{i,i} &= -r-t_{\textit{ij}} v_{(i,j)} , \end{align}
\begin{align} \rho {A}_{\textit{ij}} + G_{ijk,k} &= P_{\textit{ij}} , \end{align}
which are closed by constitutive laws.
The balance of linear momentum (6.2) is closed by the constitutive law for the stress
\begin{equation} t_{\textit{ij}} = p {\delta }_{\textit{ij}} - t^v_{\textit{ij}} , \end{equation}
where the pressure
$p$
is given in (5.19) and the viscous stress is
\begin{align} t^v_{\textit{ij}} &= \lambda (v_{l,l}) {\delta }_{\textit{ij}} + 2 \mu \left ( v_{\langle i,j\rangle } - \eta \left ( v_{\langle i,l\rangle }A_{\textit{lj}} + A_{il} v_{\langle l,j\rangle } - \frac {2}{3}\left (v_{\langle l,m\rangle }A_{\textit{lm}}\right ){\delta }_{\textit{ij}} \right ) \right ) \nonumber \\&\!\quad + \frac {2}{3} c \phi \rho T \left ( A_{\textit{ij}} - A^2_{ll} {\delta }_{\textit{ij}} + 3A^2_{\textit{ij}} - 3A^2_{ll}A_{\textit{ij}} \right )\! , \end{align}
with the viscosities
$\lambda$
and
$\mu$
are defined in (5.18), and
$c\gt 0$
is a phenomenological transport coefficient.
The balance of granular energy (6.3) is closed by the constitutive laws for the exchange of granular energy into thermal energy as
\begin{equation} r = \varphi _1 \frac {24 (1-e_p)}{\sqrt {\pi }d_p} \rho \chi \sqrt {T^3} , \varphi _1 \gt 0 . \end{equation}
The granular heat flux is
\begin{equation} Q_i = \varphi _2 \frac {4 d_p}{\sqrt {\pi }}\rho M \chi \sqrt {T^3}( {\delta }_{il} + \vartheta _2 A_{il}) (1/T)_{,l} , \varphi _2 \gt 0 , -3/2 \lt \vartheta _2 \lt 3 , \end{equation}
which yield an anisotropic granular heat conductivity.
The balance of orientation (6.4) is closed by constitutive laws for the orientation production
\begin{align} \frac {P_{\textit{ij}}}{\rho } &= \phi \left (\frac {2}{3}\left ( v_{\langle i,j\rangle } - v_{l,l}A_{\textit{ij}} - v_{l,m}A_{\textit{lm}}{\delta }_{\textit{ij}}\right ) + v_{(i,l)}A_{\textit{lj}} + A_{il}v_{(l,j)} - 2(v_{l,m}A_{\textit{lm}})A_{\textit{ij}}\right ) \nonumber \\&\!\quad - \psi d_p^{-1}\sqrt {T} A_{\textit{ij}} , \psi \gt 0 , \end{align}
and the orientation flux is
\begin{equation} G_{\textit{ijk}} = -\rho \kappa d_p \sqrt {T} A_{\textit{ij},k} , \kappa \gt 0 , \end{equation}
which is taken to be isotropic with respect to the orientation gradient.
The bulk phenomenological transport coefficients
$\{c,\eta ,\phi ,\psi ,\kappa ,\varphi _1,\varphi _2,\vartheta _2,\kappa \}$
are functions of the objective scalars
$\{\varrho ,i_1,i_2\}$
and particle properties. For spherical particles, where
$r_p=0$
and
$A_{\textit{ij}}=0$
, (6.4), (6.9) and (6.10) are not needed, and the coefficients
$\varphi _{1-2}=1$
and
$\{\vartheta _6,\kappa \}$
are not used.
The boundary conditions at impermeable wall with inward normal
$\bar {n}_i$
are
\begin{align} r^{\mathfrak{s}} &= \varphi _3 \sqrt {2/\pi }(1-e_w)\rho (1+4\chi ) \sqrt {T^3} , \varphi _3 \gt 0 , \end{align}
\begin{align} Q_l \bar {n}_l &= \varphi _4 \frac {4}{5} \rho \bar {T} \sqrt {T} - \varphi _5 \sqrt {2/\pi }(1-e_w)\rho (1+4\chi ) \sqrt {T^3} + \tilde {v}^\tau _l t^\tau _{\textit{lm}}\bar {n}_m , \varphi _{4,5}\gt 0 , \end{align}
\begin{align} t^\tau _{il} \bar {n}_l &= \varphi _6 \frac {\zeta ^2}{2\pi } \bar {T} \rho ( 1 + 4 \chi ) \sqrt {T}({\delta }_{il}+\vartheta _6 A_{il}) \tilde {v}^\tau _l , \varphi _6 \gt 0, -3/2\lt \vartheta _6\lt 3 , \end{align}
\begin{align} G_{\textit{ijl}} \bar {n}_l &= -\kappa _w d_p^{-1} \rho \sqrt {T} \, [\![ A_{\textit{ij}} ]\!]_{r_g} , \kappa _w \gt 0 , \\ & [\![ A_{\textit{ij}} ]\!]_{r_g\lt 0}=A_{\textit{ij}}-\frac {1}{3}\bar {n}_i\bar {n}_j , \notag \\ & [\![ A_{\textit{ij}} ]\!]_{r_g\gt 0}= \frac {1}{2}\left (A_{il}\bar {n}_l\bar {n}_j+\bar {n}_iA_{jl}\bar {n}_l\right )-\frac {1}{3}\left (\bar {n}_l A_{\textit{lm}}\bar {n}_m\right ){\delta }_{\textit{ij}}+\frac {1}{3}\bar {n}_i\bar {n}_j-\frac {1}{9}{\delta }_{\textit{ij}} \; \notag . \end{align}
The boundary phenomenological transport coefficients
$\{\varphi _3,\varphi _4,\varphi _5,\varphi _6,\vartheta _6,\kappa _w\}$
are functions of the objective scalars
$\{\varrho ,i_{1-6}\}$
, the particle properties and the boundary properties. For spherical particles, where
$r_p=0$
and
$A_{\textit{ij}}=0$
, (6.14) is not used and the coefficients
$\varphi _{3-6}=1$
.
7. Conclusion
In this paper, we apply non-equilibrium thermodynamics, specifically the principles of linear irreversible thermodynamics, to construct a general framework for oriented granular gases. The state variables considered in (2.1) include the density, velocity, thermodynamic temperature, granular temperature and orientation tensor. The corresponding balance laws are the conservation of mass (2.3), conservation of linear momentum (2.4), balance of thermal energy (2.7), balance of granular energy (2.8) and balance of orientation (2.11). These balance laws require closure through constitutive relations for the associated fluxes, that is, exchange of granular energy and thermal energy, thermal heat flux, granular heat flux, stress, orientation production and orientation flux.
The entropy balance (2.17) and the Gibbs equation (2.18) are used to derive the entropy production. From the entropy production expression (2.26), the thermodynamic forces (3.2) and fluxes (3.3) are identified. Taking the fluxes to be linear functions of the forces introduces the matrix of phenomenological coefficients (3.4). By applying the Onsager–Casimir reciprocal relations, we obtained the general form of the matrix of phenomenological coefficients (3.8), that must be positive semidefinite to satisfy the second law of thermodynamics. The representation of the elements of the matrix of phenomenological coefficients as functions of the objective state variables are given in Appendix A.
The boundary conditions are constructed by considering the conservation of mass (4.2), conservation of linear momentum (4.3), balance of thermal energy (4.4), balance of granular energy (4.5), balance of orientation (4.6) and balance of entropy (4.7), at the boundaries. The entropy production is used to identify the boundary thermodynamic forces (4.14) and fluxes (4.15). Again, taking linear flux-force relations and applying the Onsager–Casimir reciprocal relations lead to the boundary phenomenological coefficient matrix (4.17), which must be positive semidefinite by the second law of thermodynamics. The representation of elements of the boundary phenomenological coefficient matrix as functions of the objective state variables are given in Appendix B.
This yields a general formulation which is then reduced by neglecting what may be considered secondary effects. The reduced formulation is compared with available results from the granular kinetic theory of spherical particles and available results for oriented particles. It is shown in § 5.2 that the reduced formulation well captures these results by recovering the associated transport coefficients. A summary for the governing equations of oriented granular gas is given in § 6.
The complete set of governing equations can be applied to solve general initial-valued problems of oriented granular gases. The generalization developed here includes the additional effects raised by oriented particles. Such effects include orientational dependency of energy exchange between thermal and granular, anisotropic granular heat flux, anisotropic shear viscosity, orientational production and flux at the bulk and boundaries.
The most significant contribution of this work is the rigorous formulation of a physically admissible generalization to granular gases of oriented particles which is consistent with the second law of thermodynamics. The rigorous formulation reveals the role of the orientation in the transport coefficients and identifies coupling that may be omitted otherwise. Specifically, the identification of the required coupling of the viscous stress,
$t^v_{\textit{ij}}$
, and the orientation,
$A_{\textit{ij}}$
, through
$C^4_{\textit{ijklm}}$
in (3.8) is neglected in the literature (Berzi et al. Reference Berzi, Thai-Quang, Guo and Curtis2016; Hidalgo et al. Reference Hidalgo, Szabó, Gillemot, Börzsönyi and Weinhart2018; Nadler Reference Nadler2021; Vescovi et al. Reference Vescovi, Nadler and Berzi2024) as these generalizations suggest adding correction terms that are based on intuition rather than rigorous formulation.
To the best of our knowledge, these transport coefficients remain largely unknown. Therefore, this framework provides direction for future work to determine these coefficients either from microscopic theories or by extracting them phenomenologically from discrete microscopic simulations or experiments.
Funding
This work was supported by the Natural Sciences and Engineering Research Council of Canada through Discovery Grant RGPIN-2022-03188 (Struchtrup) and RGPIN-2018-04573 (Nadler).
Declaration of interests
The authors declare no conflict of interest.
Data availability statement
Not applicable.
Author contributions
Conceptualization, M.A., H.S. and B.N.; methodology, M.A., H.S. and B.N.; formal analysis, M.A., H.S. and B.N.; writing – original draft preparation, M.A. and B.N; writing – review and editing, M.A., H.S. and B.N.; supervision, B.N. and H.S.; funding acquisition, B.N. and H.S. All authors have read and agreed to the published version of the manuscript.
Appendix A. Representation of the bulk phenomenological coefficients
The phenomenological coefficients in (3.8) consist of tensors of order zero
$\{\mathcal{A}^1\}$
, two
$\{\mathcal{B}_{\textit{ij}}^{1-5}\}$
, four
$\{\mathcal{C}_{\textit{ijlm}}^{1-5}\}$
, five
$\{\mathcal{D}_{\textit{ijklm}}^{1,2}\}$
and six
$\{\mathcal{E}_{ijklmn}\}$
. Their representations depend on the objective scalars
$\{\rho ,\theta ,T,i_1,i_2\}$
and are constructed by combination of the tensors
$\{{\delta }_{\textit{ij}},\epsilon _{\textit{ijk}},A_{\textit{ij}}\}$
. Their tensor structures must also account for the properties, symmetry and traceless, of the associated force (3.2) and flux (3.3). Although the tensorial representations could include terms up to sixth order in the orientation tensor, e.g.
$A_{it}A_{jt}A_{ks}A_{ls}A_{mr}A_{nr}$
for
$\mathcal{E}_{ijklmn}$
, we restrict the representations to quadratic terms in orientation tensor,
$A_{\textit{ij}}$
.
To clarify the notation, we use no overhead for coefficients independent of
$A_{\textit{ij}}$
, a check accent
$\check {\square }$
for those linear in
$A_{\textit{ij}}$
, and a hat accent
$\hat {\square }$
for those quadratic in
$A_{\textit{ij}}$
.
A.1. Tensor of order zero
The zero-order tensors are scalars that can depend only on objective scalar state variables (3.5), therefore, no tensorial representation is required.
A.2. Tensor of order two
By (3.8), the second-order tensors
$\mathcal{B}^{1,2}_{\textit{ij}}$
map a second-order tensor to a scalar and
$\mathcal{B}^{3,4,5}_{\textit{ij}}$
map a vector to a vector. By (3.8), it follows that
$\mathcal{B}^1_{\textit{ij}}=\mathcal{B}^1_{(ij)}$
, hence its structure is
\begin{equation} \mathcal{B}^1_{\textit{ij}} = b_1 {\delta }_{\textit{ij}} + \check {b}_1 A_{\textit{ij}} + \hat {b}_1 A^2_{(ij)} . \end{equation}
By (3.8), it follows that
$\mathcal{B}^2_{\textit{ij}}=\mathcal{B}^2_{\langle ij\rangle }$
, hence its structure is
\begin{equation} \mathcal{B}^2_{\textit{ij}} = \check {b}_2 A_{\textit{ij}} + \hat {b}_2 A^2_{\langle ij\rangle } , \end{equation}
which vanishes for
$A_{\textit{ij}}=0$
. By (3.8), it follows that there is no requirement on
$\mathcal{B}^{3,4,5}_{\textit{ij}}$
, hence their structure follow a similar form,
\begin{equation} \mathcal{B}^\varkappa _{\textit{ij}} = b_\varkappa {\delta }_{\textit{ij}} + \check {b}_\varkappa A_{\textit{ij}} + \hat {b}_\varkappa A^2_{\textit{ij}} , \varkappa = 3,4,5 . \end{equation}
These are the most general representations of second-order tensors (Wang Reference Wang1969), as higher-order polynomials in
$A_{\textit{ij}}$
are not required by the Cayley–Hamilton theorem.
A.3. Tensor of order four
By (3.8), the fourth-order tensors
$\{\mathcal{C}^{1,2}_{\textit{ijkl}}\}$
map a third-order tensor to a vector, and
$\{\mathcal{C}^{3-5}_{\textit{ijkl}}\}$
map a second-order tensor to a second-order tensor. By (3.8), it follows that the structures of
$\mathcal{C}^{1,2}_{\textit{ijkl}}=\mathcal{C}^{1,2}_{i \langle jk \rangle l}$
are
\begin{align} \mathcal{C}^\varkappa _{\textit{ijkl}} &= c_\varkappa {\delta }_{i\langle j}{\delta }_{k\rangle l} + \check {c}_{\varkappa _1} A_{i\langle j}{\delta }_{k\rangle l} + \check {c}_{\varkappa _2} {\delta }_{i\langle j}A_{k\rangle l}+\check {c}_{\varkappa _3} {\delta }_{il}A_{jk} + \hat {c}_{\varkappa _1} A_{i\langle j}A_{k\rangle l} \nonumber \\&\!\quad + \hat {c}_{1_2} A_{il}A_{jk} + \hat {c}_{\varkappa _3} {\delta }_{i\langle j} A^2_{k\rangle l}+ \hat {c}_{\varkappa _4} A^2_{\langle jk\rangle }{\delta }_{il} + \hat {c}_{\varkappa _5} A^2_{i\langle j}{\delta }_{k\rangle l} , \varkappa =1,2 . \end{align}
By (3.8), it follows that
$\mathcal{C}^3_{\textit{ijkl}}=\mathcal{C}^3_{(ij)(kl)}$
, hence its structure is
\begin{align} \mathcal{C}^3_{\textit{ijkl}} &= c_{3_1} {\delta }_{\textit{ij}}{\delta }_{kl} + c_{3_2} {\delta }_{(i(k}{\delta }_{l)j)} + \check {c}_{3_1} {\delta }_{\textit{ij}}A_{kl} + \check {c}_{3_2} A_{\textit{ij}}{\delta }_{kl} + \check {c}_{3_3} {\delta }_{(i(k}A_{l)j)} \nonumber \\&\!\quad + \hat {c}_{3_1} A_{\textit{ij}}A_{kl} + \hat {c}_{3_2} A_{(i(k}A_{l)j)} + \hat {c}_{3_3} {\delta }_{\textit{ij}} A^2_{(kl)} + \hat {c}_{3_4} {\delta }_{(i(k} A^2_{l)j)} + \hat {c}_{3_5} A^2_{(ij)} {\delta }_{kl} . \end{align}
By (3.8), it follows that
$\mathcal{C}^4_{\textit{ijkl}}=\mathcal{C}^4_{\langle ij\rangle (kl)}$
, hence its structure is
\begin{align} \mathcal{C}^4_{\textit{ijkl}} &= c_{4} {\delta }_{\langle i(k}{\delta }_{l)j\rangle } + \check {c}_{4_1} A_{\textit{ij}} {\delta }_{kl} + \check {c}_{4_2} A_{\langle i(k}{\delta }_{l)j\rangle } + \hat {c}_{4_1} A_{\textit{ij}}A_{kl} + \hat {c}_{4_2} A_{\langle i(k}A_{l)j\rangle } \nonumber \\&\!\quad + \hat {c}_{4_3} A^2_{\langle i(k} {\delta }_{l)j\rangle } + \hat {c}_{4_4} A^2_{\langle ij\rangle } {\delta }_{kl} . \end{align}
By (3.8), it follows that
$\mathcal{C}^5_{\textit{ijkl}}=\mathcal{C}^5_{\langle ij\rangle \langle kl\rangle }$
, hence its structure is
\begin{align} \mathcal{C}^5_{\textit{ijkl}} = c_{5} {\delta }_{\langle i\langle k}{\delta }_{l\rangle j\rangle } + \check {c}_{5} A_{\langle i\langle k}{\delta }_{l\rangle j\rangle } + \hat {c}_{5_1} A_{\textit{ij}}A_{kl} + \hat {c}_{5_2} A_{\langle i\langle k}A_{l\rangle j\rangle } + \hat {c}_{5_3} {\delta }_{\langle i\langle k} A^2_{l\rangle j\rangle } . \end{align}
For the most general representation of fourth-order tensors, third-order terms such as
$A^2_{\textit{ij}}A_{kl}$
and fourth-order terms such as
$A^2_{\textit{ij}}A^2_{kl}$
should also be included.
A.4. Tensor of order five
By (3.8), fifth-order tensors map a second-order tensor to a second-order tensor. Also by (3.8) it follows that
$\mathcal{D}^1_{\textit{ijklm}}=\mathcal{D}^1_{(ij)\langle kl\rangle m}$
, hence its structure is
\begin{align} \mathcal{D}^1_{\textit{ijklm}} &= d_1 {\delta }_{(i\langle k}\epsilon _{l\rangle j)m} + \check {d}_{1} A_{(i\langle k}\epsilon _{l\rangle j)m} +\hat {d}_{1_1} A^2_{(i\langle k}\epsilon _{l\rangle j)m} + \hat {d}_{1_2} A_{t\langle k}A_{l\rangle (i}\epsilon _{j)tm} \nonumber \\&\!\quad + \hat {d}_{1_3} A_{kl}A_{t(i}\epsilon _{j)tm} + \hat {d}_{1_4} A_{\textit{ij}}A_{t\langle k}\epsilon _{l\rangle tm} + \hat {d}_{1_5} A_{t(i}A_{j)\langle k}\epsilon _{l\rangle tm} + \hat {d}_{1_6} A_{t(i}\epsilon _{j)t\langle k}A_{l\rangle m} \nonumber \\&\!\quad + \hat {d}_{1_7} A_{t\langle k}\epsilon _{l\rangle t(i}A_{j)m}+ \hat {d}_{1_8} A_{(i\langle k}\epsilon _{l\rangle j)t}A_{tm} . \end{align}
By (3.8), it follows that
$\mathcal{D}^2_{\textit{ijklm}}=\mathcal{D}^2_{\langle ij\rangle \langle kl\rangle m}$
, hence its structure is
\begin{align} \mathcal{D}^2_{\textit{ijklm}} &= d_2 {\delta }_{\langle i\langle k}\epsilon _{l\rangle j\rangle m} + \check {d}_{2} A_{\langle i\langle k}\epsilon _{l\rangle j\rangle m} + \hat {d}_{2_1} A^2_{\langle i\langle k}\epsilon _{l\rangle j\rangle n} + \hat {d}_{2_2} A_{t\langle k}A_{l\rangle \langle i}\epsilon _{j\rangle tm} \nonumber \\&\!\quad + \hat {d}_{2_3} A_{kl}A_{t\langle i}\epsilon _{j\rangle tm} + \hat {d}_{2_4} A_{\textit{ij}}A_{t\langle k}\epsilon _{l\rangle tm} + \hat {d}_{2_5} A_{t\langle i}A_{j\rangle \langle k}\epsilon _{l\rangle tm} + \hat {d}_{2_6} A_{t\langle i}\epsilon _{j\rangle t\langle k}A_{l\rangle m} \nonumber \\&\!\quad + \hat {d}_{2_7} A_{t\langle k}\epsilon _{l\rangle t\langle i}A_{j\rangle m}+ \hat {d}_{2_8} A_{\langle i\langle k}\epsilon _{l\rangle j\rangle t}A_{tm} . \end{align}
For the most general representation of fifth-order tensors, third-order tensors such as
$A^2_{\textit{ij}}A_{kl}$
should also be included.
A.5. Tensor of order six
By (3.8), a sixth-order tensor maps a third-order tensor to a third-order tensor, it follows that
$\mathcal{E}_{ijklmn}=\mathcal{E}_{\langle ij\rangle k\langle lm\rangle n}$
, hence its structure is
\begin{align} \mathcal{E}_{ijklmn} &= e_{1} {\delta }_{k\langle i}{\delta }_{j\rangle \langle l}{\delta }_{m\rangle n} + e_2 {\delta }_{k\langle l}{\delta }_{m\rangle \langle i}{\delta }_{j\rangle n} + e_{3} {\delta }_{\langle i\langle l}{\delta }_{m\rangle j\rangle }{\delta }_{kn} + e_{4} \epsilon _{k\langle i\langle l}\epsilon _{m\rangle j\rangle n} \nonumber \\&\!\quad + \check {e}_{1} A_{k\langle i}{\delta }_{j\rangle \langle l}{\delta }_{m\rangle n}+\check {e}_{2} A_{k\langle l}{\delta }_{m\rangle \langle i}{\delta }_{j\rangle n} + \check {e}_{3}A_{\textit{ij}}{\delta }_{k\langle l}{\delta }_{m\rangle n}+\check {e}_{4}{\delta }_{k\langle i}A_{j\rangle \langle l}{\delta }_{m\rangle n} \nonumber \\&\!\quad + \check {e}_{5} {\delta }_{k\langle l}A_{m\rangle \langle i}{\delta }_{j\rangle n} + \check {e}_{6}A_{\textit{lm}}{\delta }_{k\langle i}{\delta }_{j\rangle n} + \check {e}_{7}{\delta }_{k\langle l}{\delta }_{m\rangle \langle i}A_{j\rangle n} + \check {e}_{8}{\delta }_{k\langle i}{\delta }_{j\rangle \langle l}A_{m\rangle n} \nonumber \\&\!\quad + \check {e}_{9} A_{\langle i\langle l}{\delta }_{m\rangle j\rangle }{\delta }_{kn} + \check {e}_{10} {\delta }_{\langle i\langle l}A_{m\rangle j\rangle } {\delta }_{kn} + \check {e}_{11} {\delta }_{\langle i\langle l}{\delta }_{m\rangle j\rangle }A_{kn} \nonumber \\&\!\quad +\hat {e}_{1} {\delta }_{k\langle i}A_{j\rangle \langle l}A_{m\rangle n}+\hat {e}_{2} {\delta }_{k\langle i}A_{j\rangle n}A_{\textit{lm}} + \hat {e}_{3}{\delta }_{k\langle i}A_{j\rangle \langle l}A_{m\rangle n}+\hat {e}_{4}{\delta }_{k\langle l}A_{m\rangle n}A_{\textit{ij}} \nonumber \\&\!\quad + \hat {e}_{5} {\delta }_{\langle i\langle l}A_{m\rangle j\rangle }A_{kn} + \hat {e}_{6}{\delta }_{k\langle l}A_{m\rangle \langle i}A_{j\rangle n} + \hat {e}_{7}A_{k\langle i}{\delta }_{j\rangle \langle l}A_{m\rangle n} + \hat {e}_{8}A_{k\langle l}{\delta }_{m\rangle \langle i}A_{j\rangle n} \nonumber \\&\!\quad + \hat {e}_{9}A_{\langle i\langle l}{\delta }_{m\rangle j\rangle }A_{kn} +\hat {e}_{10} A_{k\langle i}A_{\langle lm\rangle }{\delta }_{j\rangle n} + \hat {e}_{11} A_{k\langle l}A_{m\rangle \langle i} {\delta }_{j\rangle n} + \hat {e}_{12} A_{\textit{ij}}A_{\langle lm\rangle }{\delta }_{kn} \nonumber \\&\!\quad + \hat {e}_{13} A_{\textit{ij}}A_{k\langle l}{\delta }_{m\rangle n} + \hat {e}_{14} A_{k\langle i}A_{j\rangle \langle l}{\delta }_{m\rangle n} + \hat {e}_{15} A_{\langle i\langle l}A_{m\rangle j\rangle }{\delta }_{kn} + \hat {e}_{16} A^2_{\langle ij\rangle }{\delta }_{k\langle l}{\delta }_{m\rangle n} \nonumber \\&\!\quad + \hat {e}_{17} A^2_{k\langle i}{\delta }_{j\rangle \langle l}{\delta }_{m\rangle n} + \hat {e}_{18} {\delta }_{k\langle i}A^2_{j\rangle \langle l}{\delta }_{m\rangle n} + \hat {e}_{19} {\delta }_{k\langle l}{\delta }_{m\rangle \langle i}A^2_{j\rangle n}{\delta }_{m\rangle n} + \hat {e}_{20} A^2_{\langle lm\rangle }{\delta }_{k\langle i}{\delta }_{j\rangle n} \nonumber \\&\!\quad + \hat {e}_{21} A^2_{k\langle l}{\delta }_{m\rangle \langle i}{\delta }_{j\rangle n} + \hat {e}_{22} {\delta }_{k\langle i}{\delta }_{j\rangle \langle l}A^2_{m\rangle n} + \hat {e}_{23} {\delta }_{\langle i\langle l}{\delta }_{m\rangle j\rangle }A^2_{km} , \end{align}
which is not the most general representation as polynomial of orders three, four, five and six of
$A_{\textit{ij}}$
are available.
Appendix B. Representation of the boundary phenomenological coefficients
The phenomenological coefficients in (4.17) consist of tensors of order zero
$\{\mathcal{A}^{2-7}\}$
, one
$\{\mathcal{H}_{i}^{1-3}\}$
, two
$\{\mathcal{B}_{\textit{ij}}^{6-9}\}$
, three
$\{\mathcal{G}_{\textit{ijk}}\}$
and four
$\{\mathcal{C}_{\textit{ijkl}}^{6}\}$
. Their representations are constructed by combinations of the tensors
$\{n_i,{\delta }_{\textit{ij}},\epsilon _{\textit{ijk}},A_{\textit{ij}},\bar {A}_{\textit{ij}}\}$
and must depend on the objective scalars (4.19). Their tensorial structures must also account for the properties, symmetry and traceless, of the associated force (4.14) and flux (4.15). Although the tensorial representations could include terms up to fourth order in the orientation tensor, e.g.
$A_{it}A_{jt}A_{ks}A_{ls}$
for
$\mathcal{C}^6_{\textit{ijkl}}$
, we restrict the representations to quadratic terms in
$A_{\textit{ij}}$
.
To clarify the notation, we use no overhead for transport coefficients that are independent of the orientation tensors
$A_{\textit{ij}}$
,
$\check {\square }$
for those linear in
$A_{\textit{ij}}$
, and
$\hat {\square }$
for those quadratic in
$A_{\textit{ij}}$
.
For convenience in constructing these representations, we introduce the projection
\begin{equation} \mathbb{P}_{\textit{ij}} = {\delta }_{\textit{ij}} - \bar {n}_i \bar {n}_j , \end{equation}
which projects a vector to the plan with normal
$\bar {n}_i$
.
B.1. Tensor of order zero
The zero-order tensors are scalars that can depend only on the objective scalar state variables, therefore, no tensorial representation is required.
B.2. Tensor of order one
By (4.17), the boundary tensors of order one,
$\mathcal{H}_i^{1,2,3}$
, map a tangent vector to a scalar. To eliminate any arbitrary component that does not contribute to the representations of such mappings, they must satisfy the constraint
$\mathcal{H}_i \bar {n}_i=0$
. Hence, the structures of
$\mathcal{H}^{1,2,3}_i$
are similar and are given by
\begin{align} \mathcal{H}^\varkappa _i &= h^\varkappa _1\mathbb{P}_{il}\bar {A}_{\textit{lm}} \bar {n}_m\! + h^\varkappa _2\mathbb{P}_{il}\bar {A}_{\textit{lm}}^2 \bar {n}_m\! + \check {h}^\varkappa _1\mathbb{P}_{il}A_{\textit{lm}} \bar {n}_m\! + \check {h}^\varkappa _2\mathbb{P}_{il}A_{\textit{lm}}\bar {A}_{\textit{mn}} \bar {n}_n + \check {h}^\varkappa _3\mathbb{P}_{il}\bar {A}_{\textit{lm}}A_{\textit{mn}} \bar {n}_n \nonumber \\&\!\quad + \check {h}^\varkappa _5\mathbb{P}_{il}\bar {A}_{\textit{lm}} A_{\textit{mn}}\bar {A}_{\textit{np}} \bar {n}_p + \check {h}^\varkappa _6\mathbb{P}_{il}\bar {A}^2_{\textit{lm}}A_{\textit{mn}}\bar {n}_n + \check {h}^\varkappa _4\mathbb{P}_{il}A_{\textit{lm}}\bar {A}^2_{\textit{mn}} \bar {n}_n + \hat {h}^\varkappa _2\mathbb{P}_{il}A^2_{\textit{lm}}\bar {A}_{\textit{mn}} \bar {n}_n \nonumber \\&\!\quad + \hat {h}^\varkappa _3\mathbb{P}_{il}\bar {A}_{\textit{lm}}A^2_{\textit{mn}} \bar {n}_n + \hat {h}^\varkappa _4\mathbb{P}_{il}A^2_{\textit{lm}}\bar {A}^2_{\textit{mn}} \bar {n}_n + \hat {h}^\varkappa _1\mathbb{P}_{il}A^2_{\textit{lm}} \bar {n}_m\! + \hat {h}^\varkappa _5\mathbb{P}_{il}\bar {A}^2_{\textit{lm}}A^2_{\textit{mn}} \bar {n}_n \nonumber \\&\!\quad + \hat {h}^\varkappa _6\mathbb{P}_{il}\bar {A}_{\textit{lm}} A^2_{\textit{mn}}\bar {A}_{\textit{np}} \bar {n}_p +\hat {h}^\varkappa _7\mathbb{P}_{il}(A\bar {A})^2_{\textit{lm}} \bar {n}_m\! + \hat {h}^\varkappa _8\mathbb{P}_{il}(\bar {A}A)^2_{\textit{lm}} \bar {n}_m , \varkappa =\{1,2,3\} . \end{align}
It should be noted that these are the most general representations (Wang Reference Wang1969) of such first-order tensors, as higher-order polynomials in the orientation tensor are not required by the Cayley–Hamilton theorem.
B.3. Tensor of order two
By (4.17), the boundary tensors
$\mathcal{B}^{6,7,8}_{\textit{ij}}$
of order two map a second-order tensor to a scalar and
$\mathcal{B}^9_{\textit{ij}}$
maps a tangent vector to a tangent vector. By (4.17), it follows that
$\mathcal{B}^{6,7,8}_{\textit{ij}}=\mathcal{B}^{6,7,8}_{\langle ij\rangle }$
, hence, their structure is similar and is represented by
\begin{align} \mathcal{B}^\varkappa _{\textit{ij}} &= b_{\varkappa _1} \bar {A}_{\textit{ij}} + b_{\varkappa _2} \bar {A}^2_{\langle ij\rangle } + b_{\varkappa _3} \bar {n}_{\langle i} \bar {n}_{j\rangle } + b_{\varkappa _4} \bar {n}_{l}\bar {A}_{l\langle i} \bar {n}_{j\rangle } + b_{\varkappa _5} \bar {n}_{l}\bar {A}^2_{l\langle i} \bar {n}_{j\rangle } \nonumber \\&\!\quad + \check {b}_{\varkappa _1} A_{\textit{ij}} + \check {b}_{\varkappa _2} A_{l\langle i}\bar {A}_{j\rangle l} + \check {b}_{\varkappa _3} \bar {n}_{\langle i} A_{j\rangle l} \bar {n}_l + \check {b}_{\varkappa _4} A_{l\langle i}\bar {A}_{j\rangle m}\bar {A}_{\textit{ml}} + \check {b}_{\varkappa _5} \bar {A}_{l\langle i}\bar {A}_{j\rangle m}A_{\textit{ml}} \nonumber \\&\!\quad + \check {b}_{\varkappa _6} \bar {n}_{\langle i} A_{j\rangle l}\bar {A}_{\textit{lm}} \bar {n}_m\! + \check {b}_{\varkappa _7} \bar {n}_{\langle i} \bar {A}_{j\rangle l}A_{\textit{lm}} \bar {n}_m\! + \check {b}_{\varkappa _8} \bar {n}_{\langle i} A_{j\rangle l}\bar {A}^2_{\textit{lm}} \bar {n}_m\! + \check {b}_{\varkappa _9} \bar {n}_{\langle i} \bar {A}_{j\rangle l}\bar {A}_{\textit{lm}} A_{\textit{mn}} \bar {n}_n \nonumber \\&\!\quad + \check {b}_{\varkappa _{10}} \bar {n}_{\langle i} \bar {A}_{j\rangle l} A_{\textit{lm}} \bar {A}_{\textit{mn}} \bar {n}_n + \hat {b}_{\varkappa _1} A^2_{\langle ij\rangle } + \hat {b}_{\varkappa _2} A^2_{l\langle i}\bar {A}_{j\rangle l} + \hat {b}_{\varkappa _3} A^2_{l\langle i}\bar {A}^2_{j\rangle l} + \hat {b}_{\varkappa _4} A^2_{l\langle i}\bar {A}^2_{j\rangle l} \nonumber \\&\!\quad + \hat {b}_{\varkappa _5} A_{l\langle i}A_{j\rangle m}\bar {A}_{\textit{ml}} + \hat {b}_{\varkappa _6} A_{l\langle i}A_{j\rangle m}\bar {A}^2_{\textit{ml}} + \hat {b}_{\varkappa _7} \bar {A}_{l\langle i}\bar {A}_{j\rangle m}A^2_{\textit{ml}} + \hat {b}_{\varkappa _7} A_{l\langle i}\bar {A}_{j\rangle m}A_{\textit{mn}}\bar {A}_{nl} \nonumber \\&\!\quad + \hat {b}_{\varkappa _8} A_{l\langle i}\bar {A}_{j\rangle m}\bar {A}_{\textit{mn}}A_{nl} + \hat {b}_{\varkappa _9} \bar {n}_l A^2_{l\langle i} \bar {n}_{j\rangle } + \hat {b}_{\varkappa _{10}} \bar {n}_l \bar {A}_{\textit{lm}}A^2_{m\langle i} \bar {n}_{j\rangle } + \hat {b}_{\varkappa _{11}} \bar {n}_l A^2_{l\langle i}\bar {A}_{j\rangle m} \bar {n}_{m} \nonumber \\&\!\quad + \hat {b}_{\varkappa _{12}} \bar {n}_l \bar {A}^2_{\textit{lm}}A^2_{m\langle i} \bar {n}_{j\rangle } + \hat {b}_{\varkappa _{13}} \bar {n}_l A^2_{l\langle i}\bar {A}^2_{j\rangle m} \bar {n}_{m} , \varkappa =\{6,7,8\} . \end{align}
By (4.17), it follows that
$\mathcal{B}^9_{il}\bar {n}_l=0$
and
$\mathcal{B}^9_{li}\bar {n}_l=0$
, hence, its structure is
\begin{align} \mathcal{B}^9_{\textit{ij}} &= b_{9_1} \mathbb{P}_{il}\bar {A}_{\textit{lm}}\mathbb{P}_{\textit{mj}} + b_{9_2}\mathbb{P}_{il}\bar {A}^2_{\textit{lm}}\mathbb{P}_{\textit{mj}} + b_{9_2}\mathbb{P}_{il}\bar {A}_{\textit{lm}} \bar {n}_m \bar {n}_n \bar {A}_{\textit{np}}\mathbb{P}_{\textit{pj}}+ \check {b}_{9_1} \mathbb{P}_{il}A_{\textit{lm}}\mathbb{P}_{\textit{mj}} \nonumber \\&\!\quad + \check {b}_{9_2} \mathbb{P}_{il}A_{\textit{lm}}\bar {A}_{\textit{mn}}\mathbb{P}_{nj} + \check {b}_{9_3} \mathbb{P}_{il}A_{\textit{lm}}\bar {A}^2_{\textit{lm}}\mathbb{P}_{\textit{mj}} + \check {b}_{9_4} \mathbb{P}_{il}\bar {A}_{\textit{lm}}A_{\textit{mn}}\bar {A}_{\textit{np}}\mathbb{P}_{\textit{pj}} \nonumber \\&\!\quad + \check {b}_{9_5} \mathbb{P}_{il}A_{\textit{lm}} \bar {n}_m \bar {n}_n\bar {A}_{\textit{np}}\mathbb{P}_{\textit{pj}} + \check {b}_{9_6} \mathbb{P}_{il}\bar {A}_{\textit{lm}} \bar {n}_m \bar {n}_n A_{\textit{np}}\mathbb{P}_{\textit{pj}} + \check {b}_{9_7} \mathbb{P}_{il} A_{\textit{lm}} \bar {n}_m \bar {n}_n \bar {A}^2_{\textit{np}}\mathbb{P}_{\textit{pj}} \nonumber \\&\!\quad + \check {b}_{9_8} \mathbb{P}_{il}\bar {A}^2_{\textit{lm}} \bar {n}_m \bar {n}_n A_{\textit{np}}\mathbb{P}_{\textit{pj}} + \hat {b}_{9_1} \mathbb{P}_{il}A^2_{\textit{lm}}\mathbb{P}_{\textit{mj}} + \hat {b}_{9_2} \mathbb{P}_{il} A^2_{\textit{lm}}\bar {A}_{\textit{mn}}\mathbb{P}_{nj} + \hat {b}_{9_3} \mathbb{P}_{il}\bar {A}_{\textit{lm}}A^2_{\textit{mn}}\mathbb{P}_{nj}\nonumber \\&\!\quad + \hat {b}_{9_4} \mathbb{P}_{il}A^2_{\textit{lm}}\bar {A}^2_{\textit{mn}}\mathbb{P}_{nj} + \hat {b}_{9_5} \mathbb{P}_{il}\bar {A}^2_{\textit{lm}}A^2_{\textit{mn}}\mathbb{P}_{nj} + \hat {b}_{9_6} \mathbb{P}_{il}\bar {A}_{\textit{lm}}A^2_{\textit{mn}}\bar {A}_{\textit{np}}\mathbb{P}_{\textit{pj}} \nonumber \\&\!\quad + \hat {b}_{9_7} \mathbb{P}_{il}A^2_{\textit{lm}} \bar {n}_m \bar {n}_n\bar {A}_{\textit{np}}\mathbb{P}_{\textit{pj}} + \hat {b}_{9_8} \mathbb{P}_{il}\bar {A}_{\textit{lm}} \bar {n}_m \bar {n}_n A^2_{\textit{np}}\mathbb{P}_{\textit{pj}} + \hat {b}_{9_9} \mathbb{P}_{il}A^2_{\textit{lm}} \bar {n}_m \bar {n}_n \bar {A}^2_{\textit{np}}\mathbb{P}_{\textit{pj}} \nonumber \\&\!\quad + \hat {b}_{9_{10}} \mathbb{P}_{il}\bar {A}^2_{\textit{lm}} \bar {n}_m \bar {n}_n A^2_{\textit{np}}\mathbb{P}_{\textit{pj}} + \check {b}_{9_{11}} \mathbb{P}_{il}A_{\textit{lm}}\bar {A}_{\textit{mn}} \bar {n}_n \bar {n}_p A_{pq}\bar {A}_{qs}\mathbb{P}_{ps} \nonumber \\&\!\quad + \check {b}_{9_{12}} \mathbb{P}_{il}\bar {A}_{\textit{lm}}A_{\textit{mn}} \bar {n}_n \bar {n}_p A_{pq} \bar {A}_{qs}\mathbb{P}_{sj} + \check {b}_{9_{13}} \mathbb{P}_{il}A_{\textit{lm}}\bar {A}_{\textit{mn}} \bar {n}_n \bar {n}_p\bar {A}_{pq}A_{qs}\mathbb{P}_{sj} \nonumber \\&\!\quad + \check {b}_{9_{14}} \mathbb{P}_{il}\bar {A}_{\textit{lm}}A_{\textit{mn}} \bar {n}_n \bar {n}_p \bar {A}_{pq} A_{qs}\mathbb{P}_{sj} . \end{align}
B.4. Tensor of order three
By (4.17), the boundary third-order tensor
$\mathcal{G}_{\textit{ijk}}$
maps a second-order tensors to a tangent vector, hence, it must satisfy
$\mathcal{G}_{ljk} \bar {n}_l=0$
in addition to the property
$\mathcal{G}_{\textit{ijk}} =\mathcal{G}_{i\langle jk\rangle }$
.
Hence, its structure is
\begin{align} \mathcal{G}_{\textit{ijk}} &= g_1\mathbb{P}_{il}{\delta }_{l\langle j}\bar {n}_{k\rangle } + g_2\mathbb{P}_{il}\bar {A}_{l\langle j}\bar {n}_{k\rangle } + g_3\mathbb{P}_{il}\bar {A}^2_{l\langle j}\bar {n}_{k\rangle } + g_4\mathbb{P}_{il}\bar {A}_{\textit{lm}} \bar {n}_m \bar {n}_{\langle j}\bar {n}_{k\rangle } \nonumber \\&\!\quad + g_5\mathbb{P}_{il}\bar {A}^2_{\textit{lm}} \bar {n}_m \bar {n}_{\langle j}\bar {n}_{k\rangle } +\check {g}_1\mathbb{P}_{il} A_{l\langle j}\bar {n}_{k\rangle } + \check {g}_2\mathbb{P}_{il}A_{l\langle j}\bar {A}_{k\rangle m}\bar {n}_{m} + \check {g}_3\mathbb{P}_{il}\bar {A}_{l\langle j}A_{k\rangle m} \bar {n}_{m} \nonumber \\&\!\quad + \check {g}_4\mathbb{P}_{il}A_{l\langle j}\bar {A}^2_{k\rangle m}\bar {n}_{m} + \check {g}_5\mathbb{P}_{il}\bar {A}^2_{l\langle j}A_{k\rangle m} \bar {n}_m\! + \check {g}_6\mathbb{P}_{il}A_{\textit{lm}}\bar {A}_{m\langle j}\bar {A}_{k\rangle n}\bar {n}_n \nonumber \\&\!\quad + \check {g}_7\mathbb{P}_{il} A_{\textit{lm}} \bar {n}_m \bar {n}_{\langle j} \bar {n}_{k\rangle } + \check {g}_8\mathbb{P}_{il} A_{\textit{lm}}\bar {A}_{\textit{mn}} \bar {n}_n \bar {n}_{\langle j} \bar {n}_{k\rangle } + \check {g}_9\mathbb{P}_{il} \bar {A}_{\textit{lm}}A_{\textit{mn}} \bar {n}_n \bar {n}_{\langle j} \bar {n}_{k\rangle } \nonumber \\&\!\quad + \check {g}_{10}\mathbb{P}_{il} A_{\textit{lm}}\bar {A}^2_{\textit{mn}} \bar {n}_n \bar {n}_{\langle j} \bar {n}_{k\rangle } + \check {g}_{11}\mathbb{P}_{il} \bar {A}^2_{\textit{lm}}A_{\textit{mn}} \bar {n}_n \bar {n}_{\langle j} \bar {n}_{k\rangle } + \check {g}_{12}\mathbb{P}_{il} \bar {A}_{\textit{lm}}A_{\textit{mn}}\bar {A}_{\textit{np}} \bar {n}_p \bar {n}_{\langle j} \bar {n}_{k\rangle } \nonumber \\&\!\quad + \hat {g}_1\mathbb{P}_{il}A^2_{l\langle j} \bar {n}_{k\rangle } + \hat {g}_2\mathbb{P}_{il}A^2_{l\langle j}\bar {A}_{k\rangle m} \bar {n}_{m} + \hat {g}_3\mathbb{P}_{il}\bar {A}_{l\langle j}A^2_{k\rangle m} \bar {n}_{m} + \hat {g}_4\mathbb{P}_{il}A^2_{l\langle j}\bar {A}^2_{k\rangle m} \bar {n}_{m} \nonumber \\&\!\quad + \hat {g}_5\mathbb{P}_{il}\bar {A}^2_{l\langle j}A^2_{k\rangle m} \bar {n}_m + \hat {g}_6\mathbb{P}_{il}A^2_{\textit{lm}}\bar {A}_{m\langle j}\bar {A}_{k\rangle n} \bar {n}_n + \hat {g}_7\mathbb{P}_{il}A_{\textit{lm}}\bar {A}_{m\langle j}A_{k\rangle n}\bar {A}_{\textit{np}} \bar {n}_p \nonumber \\&\!\quad + \hat {g}_8\mathbb{P}_{il}\bar {A}_{\textit{lm}}A_{m\langle j}\bar {A}_{k\rangle n}A_{\textit{np}} \bar {n}_p + \hat {g}_9\mathbb{P}_{il}\bar {A}_{\textit{lm}}A_{m\langle j}A_{k\rangle nr}\bar {A}_{\textit{np}} \bar {n}_p \nonumber \\&\!\quad + \hat {g}_{10}\mathbb{P}_{il}A_{\textit{lm}}\bar {A}_{m\langle j}\bar {A}_{k\rangle n}A_{\textit{np}} \bar {n}_p + \hat {g}_{11}\mathbb{P}_{il} A^2_{\textit{lm}} \bar {n}_m \bar {n}_{\langle j} \bar {n}_{k\rangle } + \hat {g}_{12}\mathbb{P}_{il} A^2_{\textit{lm}}\bar {A}_{\textit{mn}} \bar {n}_n \bar {n}_{\langle j} \bar {n}_{k\rangle } \nonumber \\&\!\quad + \hat {g}_{13}\mathbb{P}_{il} A_{\textit{lm}}\bar {A}_{\textit{mn}}A_{\textit{np}} \bar {n}_p \bar {n}_{\langle j} \bar {n}_{k\rangle } + \hat {g}_{14}\mathbb{P}_{il} \bar {A}_{\textit{lm}}A^2_{\textit{mn}} \bar {n}_n \bar {n}_{\langle j} \bar {n}_{k\rangle } + \hat {g}_{15}\mathbb{P}_{il} A^2_{\textit{lm}}\bar {A}^2_{\textit{mn}} \bar {n}_n \bar {n}_{\langle j}\bar {n}_{k\rangle } \nonumber \\&\!\quad + \hat {g}_{16}\mathbb{P}_{il} A_{\textit{lm}}\bar {A}_{\textit{mn}}A_{\textit{np}}\bar {A}_{pq} \bar {n}_q \bar {n}_{\langle j}\bar {n}_{k\rangle } + \hat {g}_{17}\mathbb{P}_{il} A_{\textit{lm}}\bar {A}^2_{\textit{mn}}A_{ps} \bar {n}_n \bar {n}_{\langle i} \bar {n}_{j\rangle } \nonumber \\&\!\quad + \hat {g}_{18}\mathbb{P}_{il} \bar {A}^2_{\textit{lm}}A^2_{\textit{mn}} \bar {n}_n \bar {n}_{\langle i} \bar {n}_{j\rangle } + \hat {g}_{19}\mathbb{P}_{il} \bar {A}_{\textit{lm}}A_{\textit{mn}}\bar {A}_{\textit{np}}A_{pq} \bar {n}_q \bar {n}_{\langle i} \bar {n}_{j\rangle } \nonumber \\&\!\quad + \hat {g}_{20}\mathbb{P}_{il} \bar {A}_{\textit{lm}}A^2_{\textit{mn}}\bar {A}_{\textit{np}} \bar {n}_p \bar {n}_{\langle i} \bar {n}_{j\rangle } . \end{align}
It should be noted that this is (Wang Reference Wang1969) the most general representations of third-order tensors, as higher-order polynomials of
$A_{\textit{ij}}$
are not required by the Cayley–Hamilton theorem.
B.5. Tensor of order four
By (4.17), the boundary fourth-order tensor,
$\mathcal{C}^6_{\textit{ijkl}}$
, maps a second-order tensors into a second-order tensor and it has the property
$\mathcal{C}^6_{\textit{ijkl}}=\mathcal{C}^6_{\langle ij\rangle \langle kl\rangle }$
. Hence, its structure is
Open access
