New Foundations of Newman’s Theory for Solid Electrolytes: Thermodynamics and Transient Balances

Newman’s developments of electrochemical thermodynamics and transport phenomena are extended to account for changes in elastic or plastic deformation stress, as well as for varying pressure, temperature, and composition. The material, momentum, energy, and entropy balances from concentrated-solution transport theory are modified to include extended underpinning thermodynamic relationships. Various colligative properties describing electrolytic materials at equilibrium are identified through analysis of the Gibbs function, which also gives rise to state equations for volume, deformation strain, and entropy. A dynamical analysis on the continuum scale produces multiphysical versions of the equation of motion and the heat equation. The discussion concludes with a statement of the second law from the perspective of irreversible thermodynamics, providing a basis for the future development of transport constitutive laws for ion conductors with elastic, as well as viscous, character. © The Author(s) 2017. Published by ECS. This is an open access article distributed under the terms of the Creative Commons Attribution 4.0 License (CC BY, http://creativecommons.org/licenses/by/4.0/), which permits unrestricted reuse of the work in any medium, provided the original work is properly cited. [DOI: 10.1149/2.0611711jes] All rights reserved.

Newman's concentrated-solution theory 1 provides one of the most general continuum-scale frameworks for describing the motion of dissolved ions through liquids in response to applied composition gradients or electric fields. Some time ago Newman extended the theory to include thermoelectric effects when temperature, as well as composition, varies locally. 2 Monroe and Delacourt recently reported an extension for systems where local electroneutrality is violated, 3 which also formalized the algorithm Newman used to develop fluxexplicit transport laws for binary electrolytes in his seminal textbook 1 and extended it to multicomponent situations.
In response to a simplification made in the landmark battery model of Doyle,Fuller,and Newman,4,5 Liu and Monroe studied phenomena arising from considering the variation of an electrolyte's density with respect to its local composition. 6,7 They found that application of density's dependence on composition as a local constraint can substantially affect transport, particularly when electrolytes are non-aqueous. In simulations with multiple spatial dimensions, it appears that when a volume-explicit equation of state is adopted to describe density, model closure requires a momentum balance to be carried alongside the material balances. 6 This begs for an investigation of effects arising from stress, which have not been studied in detail within Newman's theory.
Weber and Newman 8 have analyzed the deformation stresses that result from changes in the hydration levels of polymeric fuel-cell membranes. In their study, stress was computed in a global way, being derived from a membrane's total water content. Meyers, and later, Weber and Newman, considered local mechanical forces in the form of a pressure-gradient driving force for water transport; 9,10 this pressure, however, was not connected to the nonuniform volume distribution arising from water-content gradients, which could also translate into nonuniform stress. It would be useful to handle stresses in a local way, while also addressing the coupling that potentially exists between mass and momentum transport.
This paper aims to lay out the fundamental balance equations on which concentrated-solution theory is based in a way that permits a consistent inclusion of local pressure or deformation stress. The proposed model will be used in future work to identify transport constitutive laws applicable to elastic, ionically conductive electrolytic media such as ceramic electrolytes or swollen polymer gels. * Electrochemical Society Student Member. z E-mail: priyamvada.goyal@eng.ox.ac.uk; charles.monroe@eng.ox.ac.uk

Advantages of Newman's Approach
The rigor of concentrated-solution theory owes to its grounding in classical equilibrium thermodynamics and the Onsager-Stefan-Maxwell (OSM) transport formalism, the latter of which constitutes a sort of 'irreversible thermodynamics'. OSM diffusion driving forces emerge from a Gibbs-Duhem relation among intensive-property gradients; the fluxes conjugate to the OSM forces are identified by analyzing the irreversible energy dissipation associated with dynamic processes. [11][12][13] Since Newman's theory is grounded in classical thermodynamics, it involves familiar equilibrium material properties, as well as ensuring that those properties depend on the same basis of intensive variables needed to specify an equilibrium state. More specifically, properties such as volumetric heat capacity, isothermal compressibility, and molar volume appear in Newman's theory; these can all in principle be extracted from isolated measurements of macroscopic systems at fixed composition and ambient conditions (temperature, pressure, etc.).
When applied to materials wherein properties vary locally, Newman's theory incorporates two assumptions of 'microscopic equilibrium'. First, the material properties in nonuniform dynamical systems are taken to depend pointwise on the instantaneous local values of the thermodynamic basis variables. Second, the local property dependences are assumed to identify with those observed for uniform, macroscopically equilibrated systems. These continuum-level hypotheses agree with those Onsager expressed when deriving his reciprocal relation. 13 Transport properties in Newman's theory also have explicit roots in nonequilbrium statistical mechanics. Irreversible thermodynamics as deployed by Newman connects to microscopic statistics through fluctuation theory, which is supported by a regression hypothesis developed in broad strokes by Onsager, and more practically by Casimir. 11,14 When applied to the OSM multicomponent diffusion model, fluctuation theory not only proves a symmetric reciprocal relation among the transport coefficients, but also justifies a statement that the transport coefficients depend on the same basis of thermodynamic state variables as the equilibrium material properties. 15 Wheeler and Newman have also drawn these conclusions without an overt regression hypothesis, from a nonequilibrium statistical-mechanical approach based on the fluctuation-dissipation theorem and Nyquist-Johnson linear-response theory. 16,17 A unique feature of the OSM formalism is its ability to account for 'coupled diffusion' processes, in which a driving force induces fluxes to which it is not conjugate in the function describing local energy dissipation. Diffusion in binary electrolytes typefies a coupled process, because the gradient of a given species' electrochemical potential can generally drive the flux of a different species. (A lucid example is provided by electro-osmotic drag in low-temperature fuel-cell membranes: an electric field drives migration of mobile ions relative to the membrane backbone; these ions exert a drag force on dissolved water that causes it to move relative to the membrane as well; thus electric fields apparently drive water flux, despite the fact that water molecules carry no net charge. 18 ) This physical content distinguishes OSM theory from the more common Nernst-Planck transport formalism, which ignores solute/solute interactions entirely. Newman's theory is more 'complete' than Nernst-Planck theory in the sense that it incorporates one transport coefficient for every pair of chemically distinct diffusing species, in line with the principles of irreversible thermodynamics laid out by Onsager. 19,20 Coupled diffusion is illustrated more starkly in systems where multiple types of transport occur simultaneously, whose analysis is entirely outside the scope of Nernst-Planck theory. An example is provided by materials in which temperature gradients induce charge flux or potential gradients induce heat flux, which have been of interest to electrochemists for decades. 21 Newman has extended concentratedsolution theory to rationalize this 'Thompson-Peltier effect' by incorporating terms associated with Soret and Dufour fluxes, in accord with the principles of irreversible thermodynamics. 2 A theory that incorporates the phenomenon of pressure diffusion with similar depth and rigor would be useful.
Modern electrochemical engineering applications often require modeling of 'multiphysical' transport processes, where the assumption of isothermal, isobaric conditions needs to be relaxed. Multiphysics leads to another, subtler type of coupling: fluxes across a material can be affected by the local variation of material properties with respect to the distributions of thermodynamic basis variables. This coupling is consequently driven primarily by the general balance equations, rather than being explicitly seen in the transport constitutive laws. For example, when a current passes through a liquid electrolyte, Joule heating (or 'I 2 R' heating) can cause substantial local enthalpy generation, 22 which, when distributed according to the local thermalenergy balance, induces temperature gradients. If, say, the electrolyte diffusivity depends strongly on temperature, then the local change in diffusivity with temperature will augment composition-gradientdriven material fluxes in proportion to the temperature gradienta thermophoretic force that is inherently macroscopic, and therefore distinct from the Soret effect. Since transport properties often depend strongly on temperature and composition, 23 this 'multiphysical coupling', although more obscure in its origins, can sometimes outweigh the diffusional coupling built into the OSM laws.
The description of solid electrolytes poses a substantial challenge because in many cases such materials have elastic character, i.e., many solids can sustain spatial deformations at equilibrium. Thus the classical-thermodynamic foundations of Newman's theory must be rebuilt, because the statement of the first law from which all the thermodynamic properties derive must account for the presence of deformation stress, in addition to pressure.
We set out here to develop elementary balance equations that can be used to produce a concentrated-solution theory applicable to both liquid and solid electrolytes. The system of balances allows for multiphysical couplings between electrochemical and mechanical phenomena, and reveals the bases of fluxes and driving forces that should be involved in constitutive laws that account for diffusional coupling.
Deformation stresses are shown to propagate into the Gibbs-Duhem relation among intensive properties, and consequently to alter the familiar definition of the thermodynamic force that drives mass diffusion. Moreover, since the continuity of total energy is used to derive both the heat equation and the energy dissipation (local entropy generation), these macroscopic balances are both shown to require additional terms when describing solids. Despite the rather esoteric fundamental nature of the present enquiry, this effort is probably worthwhile, because it will permit the later formulation of thermodynamically complete descriptions of ion transport in ionomer gels and ceramic electrolytes, both of which are important classes of electrochemical materials.

Thermodynamics
Mechanical principles.-A basic axiom of classical thermodynamics holds that every equilibrium colligative property of a homogeneous phase containing n chemically distinct species at ambient temperature T and external pressure p ultimately derives from a single scalar function, called the Gibbs function G, which measures the total energy available to be extracted from the phase via reversible processes. To account for deformation equilbria in elastic solids, this elementary framework must be extended by replacing the scalar external pressure p with the more general Cauchy stress tensor σ, which quantifies how any plane's orientation at a given location maps into a force vector acting on that plane. More precisely, the traction force per unit area t exerted on a plane whose outward orientation is designated by the unit normal vector e n relates to the Cauchy stress through the projection t = − e n · σ.
(It should be noted that a number of standard operations will be used here to combine tensor or vector components, rearrange a tensor's components, or map tensor or vector quantities into scalars. Tensor operations in the sequel include the scalar product, transpose, inverse, vertical scalar product, trace, adjugate trace, and determinant. These are defined in the Appendix, which also summarizes some identities and useful properties. The appendix also addresses tensor transformations, such as those that describe rigid rotations of the coordinates under a vector or tensor's components, as well as outlining the spectral theorem, which is needed to replace tensor vertical products with products of tensor invariants.) Classical thermodynamics is generally applied to simple fluids or, more generally, to systems in which only pure compression or dilation needs to be considered when computing mechanical contributions to the Gibbs function. In this mode the change in G arising from mechanical phenomena is familiarly written as dG mech = V dp, where V is the system volume. When there is deformation stress as well as compressive stress, one has instead that where V is a tensor with units of volume that describes both the system volume and the extent to which that volume is deformed.
Being conjugate to intensive property σ, V must be extensive. Some properties of σ are required by laws of mechanics; it is worth using those constraints to put Equation 1 into a simpler form. Stresses fall into two broad classes: 'spherical' and 'deviatoric'. Spherical stresses are uniformly compressive or dilatational, manifesting as purely normal tractions whose magnitudes are unaffected by the orientations of the surfaces on which they act. Deviatoric stresses, by contrast, embody tractions with no parts that remain constant with respect to the orientations of their planes of action. Equation A10 can be used to prove that any spherical tensor is necessarily a scalar multiple of the identity tensor I . Thus pressure p is identified by decomposing the Cauchy stress into its unique spherical and deviatoric parts, Thermodynamic pressure identifies with the spherical stress, related to the first invariant of the Cauchy stress through p = 1 3 tr( σ). The deviatoric nature of the deformation stress, τ, is expressed by the property tr( τ) = 0.
The volume tensor V likewise decomposes uniquely into spherical and deviatoric parts, where V = 1 3 tr( V ) is the familiar system volume; the deformationvolume tensor has the property tr( V ) = 0 . Insertion of decompositions 2 and 3 into Equation 1 gives [4] because the vertical scalar product between spherical and deviatoric tensors generally vanishes (cf. Equation A9). Thus one recovers an expression for the change in the mechanical Gibbs function that reduces to the familiar form in the absence of deformation stress. Additional simplifications can be made to the right side of Equation 4 by considering some consequences of tensor symmetry. On physical grounds one expects that the content of a Cauchy stress should not change if the coordinate system under it is subjected to an arbitrary rigid motion. This invariance principle is equivalent to a statement that stress is symmetric, σ = σ T . 24 Stress symmetry implies that the alternative projection t = − σ · e n , as well as the original form t = − e n · σ, produces the same traction force on the surface normal to e n . Through Equation 2 it also implies symmetry of the deformation stress, so that τ = τ T .
Since the vertical scalar product between an antisymmetric tensor and a symmetric tensor vanishes, it immediately follows from the symmetry of σ that any antisymmetric part of the volume tensor will not contribute to the Gibbs energy. The requirement that all colligative equilibrium properties must derive from the Gibbs function therefore amounts to a stipulation that V = V T , which implies symmetry of the deformation-volume tensor V .
Properties of the deviators τ and V allow further simplification of Equation 4. Given a deformation stress expressed relative to a particular right-handed Cartesian coordinate system, the spectral theorem (cf. the Appendix) says that one can always rigidly rotate the coordinate frame to make the tractions represented by the stress transform into three mutually perpendicular, purely normal forces. That is, for any τ there exists a principal representation [5] in which the diagonal entries {τ 1 , τ 2 , τ 3 } of d τ establish the magnitudes of the principal deformation stresses, and the rotation Q τ transforms the original coordinate basis under τ into the principal directions { e τ 1 , e τ 2 , e τ 3 } of τ. Owing to its symmetry, the basis of the deformation-volume tensor V can similarly be rotated by some Q V to produce a principal representation that involves only the principal deformation volumes Substitution of the principal representations of d τ and V allows one to demonstrate that This quantity achieves its maximum magnitude when Q V = Q dτ , that is, when the principal axes of deformation volume and deformation stress coincide. Physically, these must coincide, or else either V or τ would have some content that does not contribute to the Gibbs function, contradicting the axiom that G accounts for all of the available energy.
The physical notion that the principal axes of deformation stress and deformation volume coincide is equivalent to the mathematical commutativity relation V · τ = τ · V . Supposing that the constitutive relation between stress and strain is such that this commutativity relation holds, it follows that [7] The principal values on the right are not independent, however, because both tensors involved are traceless. Equation A12 therefore implies that τ 3 = −(τ 1 + τ 2 ) and V 3 = −(V 1 + V 2 ); consequently From this, one can identify extensive properties that associate strains with each independent principal deformation stress: so that a simple expression involving only scalar quantities, encapsulates all the possible contributions that mechanical forces can make to changes in the Gibbs function. Finally, it is useful to specify how the deformation-volume tensor V relates to the traceless Eulerian strain deviator familiar from the continuum mechanics of solids. A suitable connection is provided by the definition making the strain deviator dimensionless. This relationship is consistent with the use of strain in Coussy's comprehensive development of thermodynamics for poroelastic systems, 25 but operates in an Eulerian, rather than a Lagrangian, frame. Working backward to reinsert the tensors on which the principal stresses and deformation volumes are based, Equation 10 can be written equivalently in terms of the strain deviator as dG mech = V dp + V : d τ, [12] a relationship that can also be obtained by direct substitution of Equation 11 into Equation 4. Since is a simple scalar multiple of V , the strain deviator can be proven to share many of the fundamental properties possessed by V . In particular, = T and · σ = σ · .
Basis of the Gibbs function.-After adopting the mechanical principles developed above, one can establish a general Gibbs function to describe a homogeneous multicomponent phase. Assume that each chemically distinct constituent i has a given total molar content within the phase, which is given the symbol n i ; in an n-constituent phase, all the molar contents can be denumerated as an n-tuple {n i } n , which we take to be an ordered list for convenience. Further we assume that through the application of appropriate experimental controls, T , p, τ 1 , τ 2 , and each of the {n i } n can be varied in isolation while holding all the other quantities constant. Mathematically speaking, the domain T , p, τ 1 , τ 2 , {n i } n comprises a basis of independent scalar variables on which the available energy depends. Each of the members of this basis set of elementary state variables is real; all but τ 1 and τ 2 are strictly non-negative.
The free energy of a homogeneous, equilibrated phase is expressed as a continuous, differentiable scalar function G(T, p, τ 1 , τ 2 , {n i } n ). The Gibbs function is single-valued, but is not an injective mapping because multiple sets of the basis variables can in principle correspond to the same available energy.
Rather than treating individual species in isolation, it is convenient to define a set of alternative variables to describe composition, reserving the single quantity Journal of The Electrochemical Society, 164 (11) E3647-E3660 (2017) to express the total molar content of a phase. With n T in hand, one can define particle fractions [14] which are also non-negative, to describe local composition. Each y i is individually bounded between zero and unity, and the whole n-tuple {y i } n satisfies i y i = 1. [15] This constraint can be thought of as a 'phase rule': in an isolated phase, one of the composition descriptors necessarily depends on the others. In other words, specifying the (n −1)-tuple {y i } n−1 is sufficient to establish the composition of an n-ary phase unambiguously. Hereafter, the term 'composition' should be taken to mean a specific set of values for the basis variables {y i } n−1 . A material with n = 1 is called 'pure'. In pure materials G(T, p, τ 1 , τ 2 , n T ) and particle fractions are unnecessary. Standard theories of thermomechanics tend to be restricted to pure materials.
Extensivity.-With definitions 13 and 14 and property 15, one can show that the transformation from the composition basis {n i } n to the basis n T , {y i } n−1 is bijective (one-to-one and onto). Thus one can equivalently write the functionality of G as The most elementary property of G is that it is an extensive state function. That is to say, given any positive scalar λ, [16] Adopting n T as the only quantity that rescales the magnitude of G in this way amounts to a physical postulate: the magnitude of the available energy within a continuous material equilibrated at fixed temperature, pressure, deformation stress, and composition varies in direct proportion to the total number of particles the material contains. The term 'extensive' can be applied to any other state variable or state function that has the property stated in Equation 16. For example, Equation 14 can be used to show that each of the state variables {n i } n is extensive, because n i (T, p, {y i } n−1 , n T ) = y i n T and λ · (y i n T ) = y i · (λn T ). In general, the sum of extensive properties proves to be extensive, so it follows through Equation 13 that n T itself is extensive.
Conversely, a property is 'intensive' if its magnitude is insensitive to changes in n T at fixed T , p, and composition. For example, the molar Gibbs energy, G, is defined as G = G/n T . By letting λ = 1/n T in Equation 16, one sees that G is an intensive state function, because its domain does not need to include n T ; specifically, Euler's theorem.-Volume V and entropy S are important characteristics of a phase, which are identified as the conjugate variables of p and T , respectively. In this context 'conjugate' means they derive from G through the partial derivatives The new extensive properties D 1 and D 2 associated with volume deformation can be seen to come from G through partial derivatives of the form By a similar process one can define electrochemical potentials μ i as conjugate variables associated with the n i , Although electrochemical potentials are given a special symbol, any property obtained by applying the thermodynamic derivative on the right of Equation 19 to an extensive property is called a 'partial molar property' of species i. Partial molar properties are notated by putting an overbar on the property from which they derive and appending a subscript to indicate the species involved in the derivative. For example, one can equivalently write μ i = G i . In a phase at fixed temperature, pressure, deformation stress, and composition, the extensivity relation expressed by Equation 16 implies that [20] which is known as Euler's extensivity theorem for the grand canonical ensemble. The extensivity property from Equation 16 can be used here to prove that the electrochemical potential μ i is intensive.
The first law.-Given the stated domain of the Gibbs function, the first law of thermodynamics can be written through definitions 17 through 19 as where the second form follows from Equations 10 and 12. By inserting and applying Equation 15, one can show with the product rule of differentiation that providing an intensive form of the first law. Since it was found that G(T, p, τ 1 , τ 2 , {y} n−1 ), Equation 22 implies Equations 22 and 23 can be useful for analyzing local energy balances. Equation 23 in particular was employed by Monroe, Wheeler, and Newman in their derivation of the symmetric reciprocal relation for isothermal, isobaric OSM diffusion. 15 Volumetric equation of state.-Volume is the most fundamental descriptor of a material's kinematic state. In multicomponent materials, one must consider that the equilibrium volume can depend on composition. To quantify the contributions that individual species make to overall material volume, one can introduce partial molar volumes The second equality here derives from a Maxwell relation, which may be easier to see if one recalls that μ i = G i . Because the molar contents are naturally independent of T , p, and τ i , differentiation of Equation 20 with respect to p at fixed T and deformation stress shows through Equations 17 and 24 that Equation 16 can be used to show that V is extensive; it follows that the properties {V i } n are intensive. Also, the extensivity of V proves that the molar volume V = V /n T is intensive. There may be other physical constraints on the functionality of V , which can be seen as fundamental constitutive relationships. An important constraint applies to the mechanical part of G as a consequence of the physical definitions of spherical and deformation E3651 stresses. Discussions of these stresses often conflate a dynamical perspective with a kinematic perspective: deformation stress is typically defined heuristically as 'stress that leads to shape change without volume change', whereas spherical stress is 'stress that leads to volume change without shape change'. These concepts are intuitive, but practically imply thermodynamic constraints. To say that τ i leads to shape change without volume change means that that is, the molar volume of a material cannot depend on deformation stress. Through the Maxwell relation it appears to follow that the deformation volumes D 1 and D 2 cannot depend on pressure p.
The total molar concentration c T is defined as the inverse of molar volume, c T = 1/V , and has the same functionality, i.e., c T (T, p, {y i } n−1 ). Thus, in elastically deformable materials, the volume-explicit equation of state needs no modification from the form used for simple fluids. Through manipulation of Equation 25 with the definitions above, any phase with properties embedded in G can be seen to obey a volumetric equation of state which provides a fundamental kinematic relation to stand alongside the phase rule provided by Equation 15. Every quantity that appears in Equation 28 depends at most on T , p, and {y i } n−1 .
If desired, one can define the molar concentration of species i at this stage, through c i = c T y i . One can also identify the volume fraction of i, which relates to its concentration and partial molar volume through φ i = V i c i . Equation 28 implies that i φ i = 1. Both molar concentrations and species volume fractions are necessarily independent of the deformation stress.
Since volume is independent of deformation stress, all differential changes in V relate to changes in T , p, and {n i } n through This equation introduces the volumetric coefficient of thermal expansion, α V , and the bulk modulus (inverse isothermal compressibility) K , two familiar material properties generally defined as [30] being intensive, α V and K are functions of temperature, pressure, and composition.
Entropic equation of state.-From the thermodynamic standpoint, a phase's entropy S can be treated similarly to its volume. Individual species' contributions to S are quantified by defining partial molar entropies S i ; these relate to electrochemical potentials through the Maxwell relation Unlike V i , S i may vary with the deformation stress.
With the domain of the Gibbs function specified as described in the paragraphs preceding Equation 13, the definition of entropy through a temperature derivative in Equation 17 shows in combination with Equations 20 and 31 that [32] Thus entropy can be allocated to individual species in the same way as the Gibbs energy and the volume.
For a system based in the Gibbs function, entropy can be obtained through temperature derivatives of G, and must consequently depend on the same basis of state variables: S(T, p, τ 1 , τ 2 , {n i } n ). The total differential of entropy can therefore be written as This includes the familiar constant-pressure heat capacity, defined as but also introduces two new thermodynamic material properties. These 'volumetric coefficients of thermal deformation' are defined by analogy to the volumetric coefficient of thermal expansion: The material properties α V , α D 1 , and α D 2 defined by Equations 30 and 35 arise in Equation 33 through which are Maxwell relations that connect the strain state to entropy. Equation 33 differs from the other thermodynamic relationships in that the terms associated with mechanical forces do not involve products of the principal deformation stress and principal deformation volume. To allow the concise formulation of dynamical balances later on, it will be convenient to embody the information contained in properties α D 1 and α D 2 into a tensor quantity, which we will call the thermal deformation tensor A D .
Define the symmetric tensor A D such that it has zero trace and commutes with the deformation-strain tensor, i.e., let [37] These characteristics respectively ensure that A D is diagonalizable, has at most two independent principal values, and is diagonalized by the same rotation tensor that diagonalizes . Since this rotation also diagonalizes τ (cf. the 'Mechanical principles' item in the Thermodynamics section), it follows that where A 1 and A 2 are the independent principal values of A D in the corresponding principal directions of τ. The entries of A D are fully specified by defining [39] In short, given the principal axes of deformation strain and the coefficients of thermal deformation, the thermal deformation tensor is formed by the similarity transformation Journal of The Electrochemical Society, 164 (11) E3647-E3660 (2017) The rotations here can also be extracted from deformation stress, Cauchy stress, etc., because the principal axes of these quantities coincide. For example, Q = Q τ .
Deformational equations of state.-Materials that can be deformed elastically require two additional kinematic state equations, which describe how the independent deformation volumes depend on intensive state variables. The discussion of the volumetric equation of state showed that typical physical definitions of compressive stress require the deformation volumes D 1 and D 2 to be independent of pressure, i.e., D k (T, τ 1 , τ 2 , {n i } n ). In light of this dependence the total differential of deformation strain D k becomes These equations introduce four material properties, given the symbols G i j because they are similar to shear moduli. Observe that the necessary pressure independence of D i requires that the quotient V /G i j is independent of pressure. Further note that there are only three independent moduli G i j because a Maxwell relation based on Equation 18, requires that G 12 = G 21 . Isotropic materials are distinguished in part by the property that 1/G 12 = 0. One can use the extensivity argument given previously to derive an Euler equation [44] Thus the deformation borne by the material as a whole can be partitioned among the individual molecular species that comprise it. This involves partial molar deformation volumes, which, like D k , must be independent of pressure. Last it should be emphasized that writing Equation 41 for both independent principal directions provides a very general equation of state for material elasticity. To see this more concretely, consider a situation where temperature and composition are constant. If the properties G i j do not vary with respect to deformation stress, and if 1/G 12 = 1/G 21 = 0 and G 11 = G 22 , then Equations 41 express Hooke's law of linear, isotropic elasticity, albeit in a differential form involving tensor invariants.

Gibbs-Duhem relations.-A canonical Gibbs-Duhem equation
arises directly from the Gibbs function, by combining the first law from Equation 21 with the total differential of Equation 20 to get i n i dμ i = −SdT + V dp + V : d τ. [45] This will be used in the Momentum balances section, because it shows how thermodynamic-property gradients associated with diffusion result in mechanical forces. (One phenomenon illustrating this connection is osmotic pressure, which arises from chemical-potential differences induced by variations in composition.) Insertion of Equation 25 into Equation 29 reveals that partial molar volumes satisfy a Gibbs-Duhem equation showing that partial molar volumes have no natural dependence on deformation stress. This result can be applied in various ways: the are both useful. Equations 25 and 46 can also be used in combination with Equation 28 to show that which expresses how changes in temperature, pressure, and intensive composition variables are linked. Last, Equations 32, 33, 38, and 40 can be combined to get a Gibbs-Duhem equation which derives from the system entropy.
Internal energy and enthalpy.-From the standpoint of irreversible thermodynamics, the constitutive laws for transport derive from inspection of the local entropy generation within a nonequilibrated material. In developing this balance one must consider the distribution of total energy, which includes additional information beyond that described by the available energy function G. Two quantities that are useful in the Energy balances section below are the enthalpy, H , which relates to G through the Legendre transformation and the internal energy U , defined as [51] Terms associated with deformation have been included in the Legendre transformation from G to U . Taking total derivatives of both sides of Equation 51 and inserting Equation 21 shows that which provides a statement of the first law that associates an extensive quantifier (entropy, volume, deformation volume, and species molar content) with every distinct mode of energy exchange (thermal, compressive, deformational, and chemical, respectively). Equation  52 can be accounted in terms of accumulation, flux, and generation; the dynamics of how extensive properties are distributed is readily described by balance equations. In locally inhomogeneous systems, extensive-property balances are stated as differential equations that govern the densities of the properties, i.e., their amounts per unit of system volume.

Specific internal energy, volume, and entropy.-Unlike the intensive variables on which G depends, the basis variables in
Energy is carried by massive particles that occupy volume, rather than being a property of the volume itself, so it is convenient to refer extensive-property densities to mass. The total mass of a system is designated by the symbol m T ; mass is commonly understood to be independent of temperature, pressure, and deformation stress. Total mass does depend on system composition, however, relating to species molar contents through where m i is the molar mass of species i. (Formally m i is a partial molar property, explaining its notation here.) The molecular hypothesis essentially states that species masses relate uniquely to their molar contents, through m i = m i n i , and that, for any given species i, m i is a constant independent of temperature, stress, and composition. Both total mass and the individual species masses are extensive.
With the total mass in hand, one can define the mass density ρ as an intensive property. Owing to the property dependence of c T established by the volumetric equation of state and the constancy of molar masses, the density proves to be determined by only temperature, pressure, and composition: ρ = ρ(T, p, {y i } n−1 ). Moreover, the quotient ρ/c T , which identifies with the total molar mass m, ρ depends on composition alone. Any extensive property can be expressed per unit of system mass to produce an intensive property. Quantities normalized in this way are referred to as 'specific' properties, and notated with a carat. For instance, specific internal energy can be written in various ways: and specific volume equals inverse mass density,V = 1/ρ. Note that the internal energy per unit volume can be written either as ρÛ or c T U . A form of the first law involving specific internal energy can be derived as follows. Use Equations 20 and 51 to express the total internal energy in terms of the state variables. Then divide this relation by m T and apply Equation 56 to get an expression forÛ in terms of fundamental intensive variables and the specific versions of their extensive conjugates. Take a total differential, then eliminate the differential intensive properties using Gibbs-Duhem Equation 45 divided by m T . Finally, insertV = 1/ρ andn i = c T y i /ρ, then use Equations 10 and 12 to show that [57] This result is needed when deriving the dynamical form of the second law in the Energy balances section. Because the equation containsŜ, a state variable that is hard to control, some further analysis is still helpful.
In particular, it is useful to have a form of Equation 57 wherein the entropy dependence is eliminated in favor of dependences on composition, temperature, and stress, which can be produced as follows. Put Equation 32 in terms of specific quantities by dividing both sides by m T ; take a total differential of this result and use Equation 49 to eliminate partial-molar-entropy changes. Thus which incorporates the thermal deformation tensor defined by Equation 40. Equations 57 and 58 combine to produce an expression for internal-energy change that does not explicitly include the entropy: This final result of thermodynamic analysis is key to the derivation of the heat equation in the Energy balances section below.

Material Balances
The macroscopic description of local dynamics within an n-ary phase begins with statements of local balance equations that express how the local time variations of properties connect to their flow. These conservation laws introduce the notions that systems change transiently and that properties differ from place to place. The most elementary expressions of continuity are universal, in the sense that the governing equations hold regardless of the particular material in question. Dynamical continuity equations complement the thermodynamic information in G and are of equivalent fundamental importance.
Material is the most basic conserved quantity. In multicomponent transport systems one must consider that each individual species within a phase obeys a continuity equation, as well as the phase as a whole. This section shows how the equations that locally express the continuity of volume, charge, and mass can be derived by combining species material balances with thermodynamic and electrical state equations. In many cases it is convenient to replace one or more of the n species balances in an n-ary system with one of these broader laws of continuity.

Molar species continuity.-The most basic governing equations
Newman uses to express the local balances of material can be written as where N i is a spatial vector representing the total molar flux of species i relative to an inertial reference frame. a (Recall that c i = c T y i .) Each constituent i within an n-ary phase obeys a local material balance in the form of Equation 60. The set of these continuity laws and the set of molar fluxes are both n-tuples, e.g., { N i } n . It is also sometimes convenient to speak in terms of velocities, rather than fluxes. These can be extracted from molar fluxes through which introduces a species velocity v i (i.e., a weighted average of individual particle velocities over all particles of type i in a volume element). The set of species velocities is also an n-tuple, { v i } n .
Volume continuity.-Alongside c T -a thermodynamic variable associated with the static volume a material occupies -it is useful to have a dynamic variable that describes how material carries volume with it as it moves. The volume-average velocity, defined as provides an elementary kinematic quantity useful to establish a local description of how a material moves in total relative to an inertial frame of reference. A macroscopic volume balance can be derived by taking the partial time derivative of Equation 48, then using the material balances from Equation 60 to eliminate terms with ∂(c T y i )/∂t. Application of the vector identity ∇ · (s v) = s ∇ · v + v · ∇s then yields The strategy used to derive this result was first presented by Newman and Chapman, whose analysis of restricted diffusion produced a volume balance for isothermal, isobaric binary electrolytes. 26 Note that one can reduce the range of the summation on the right of Equation 63 by applying the spatial identity given in Equation 47.
It is often convenient to replace one of the n local material balances with volume continuity, particularly when modeling systems where T and p are both constant. In the isothermal, isobaric case, Equation 63 provides a relationship that governs the spatial distribution of v without involving partial time derivatives. 26 Also, in isothermal, isobaric systems where the material density varies linearly with the molar concentrations of its constituents (typically a good assumption), Equation 63 simplifies dramatically, reducing to ∇ · v = 0. 6 Charge continuity.-Particle numbers can be connected to total charge, as well as total volume. Charge is an extensive quantity; the excess charge density where z i stands for the equivalent charge carried by species i and F is Faraday's constant, provides an intensive description of the static electrical state. The current density i provides a dynamic variable describing charge flow, connected to molar fluxes through Faraday's law: Multiplying Equation 60 by Fz i , then summing over all species and incorporating Equations 64 and 65, yields an equation that describes the local balance of charge. Thus, if definitions 64 and 65 are adopted, and if at least one species in a system carries a nonzero equivalent charge, one of the material balances can be replaced with which expresses charge continuity. It is common when analyzing electrochemical systems to simplify the set of material balances by replacing one of them with the charge balance. This is particularly helpful in systems subject to the condition of local electroneutrality, in which charge continuity takes the form ∇ · i = 0, allowing the instantaneous distribution of current density to be determined without considering time derivatives.
Mass continuity.-One quantity that has not yet been discussed is the species mass density, defined as [67] This can also be used to identify the mass fraction of species i, a composition variable that, like y i , is constrained such that i ω i =

Equation 55 can be inserted into Equation 68
to show that species mass fractions are similar to particle fractions in that they do not vary with temperature or stress. The total mass flux provides a dynamic variable that can be associated with convective flow of material, but it is more common to employ the mass-average velocity v. This velocity is embedded in the mass flux ρ v and relates to species fluxes as [69] Here the last equalities follow from the species velocity defined in Equation 61; ρ i v i represents the mass flux associated with the motion of species i.

Multiplication of both sides of Equation 60 by m i illustrates the continuity of individual species masses,
and summing this result across all n species demonstrates overall mass continuity, [71] Thus the elementary 'continuity equation' from fluid mechanics derives from the n material balances, and necessarily involves the massaverage velocity v.

Change of composition and flux bases.-It
should be emphasized that the n species balances expressed in the form of Equation 60 involve only n − 1 primary scalar composition variables -and that specifying the composition {y i } n−1 suffices to determine uniquely any choice of intensive composition variables (mass fractions, volume fractions, species molalities, species molarities, species mass densities, etc.) at fixed T and p.
In combination with the particle-fraction sum from Equation 15, the discussions of the volumetric equation of state, specific volume, and charge continuity above have put forward definitions that establish c T , ρ, and ρ e as functions of T , p, and {y i } n−1 . In the analysis of some transport problems, it may be convenient to invert the state equations for c T , ρ, or ρ e to exchange one of the concentrations that makes up the composition basis for the molar, mass, or charge density. For example, Newman used this strategy to solve the Poisson-Nernst-Planck equations in his analysis of current flow in polarized diffuse double layers, replacing one of the species concentrations with ρ e . 28 Similarly, there are only n primary flux variables in the n species balances, Equations 60. In some cases it is convenient to use the current density i or a bulk velocity, say, v , in place of one of the molar fluxes N i .

Convective velocities and excess fluxes.-Continuity Equation 71
is often written in a Lagrangian frame, by incorporating a convective derivative whence vector identities can be used to show that [73] The definition in Equation 72 implicitly assumes that the mass-average velocity v is the correct choice of reference velocity for convection. In fact a variety of reference velocities is available, and it is worth generally developing the many ways one can parse convective contributions to the molar flux N i . Only the part of species motion that is in excess of some 'bulk flow velocity' v ψ should be associated with dissipative modes of mass transport (diffusion, migration, etc.). One can define the excess molar flux of i relative to v ψ as to produce a set of fluxes that isolate non-convective phenomena and are therefore more natural for the formulation of transport constitutive laws.
There is some ambiguity in the definition of the 'bulk flow velocity'. The convention used when defining substantial derivatives suggests v, but we have already seen another viable candidate -the volume-average velocity v -and there are many more, including familiar examples like the mole-average velocity, 31 and less familiar ones like the enthalpy-average velocity, which Hirschfelder et al. suggest might be a natural reference for heat convection. 29 Every E3655 reference velocity is defined through a sum of the form if one calls v ψ the 'ψ-average velocity', then one understands ψ i to be an n-tuple of variables defined in terms of masses or equilibrium properties that follow from the Gibbs function; these variables should also be 'fractional', in the sense that they are constrained by a relation For instance, if one identifies the fractional variable as the mass fraction, ω i → ψ i , then Equation 75 defines the mass-average velocity v; if it is identified as the species volume fraction, φ i → ψ i , the equation defines the volume-average velocity v ; one can also identify the mole-average velocity v * by making the substitution y i → ψ i , or the enthalpy-average velocity v h by c T H i y i /ρĤ → ψ i . The choice of reference velocity generally constrains the n-tuple of excess fluxes relative to that velocity. Equations 74 and 75 can be combined to yield a kinematic relation showing that only n − 1 of the n excess fluxes relative to v ψ are independently variable. In some sense, the particular choice of convective velocity is immaterial, since the difference between any two reference velocities can be determined from either independent set of excess fluxes. For reference velocities v ψ and v ψ based on fractions ψ i and ψ i , the relations follow directly from fraction-sum 76 and kinematic relation 77. One can use this difference to show that which is an invertible, homogeneous transformation that can be used to convert any set of independent excess fluxes into any other. Monroe and Newman have shown that linear transformations of this type preserve the symmetry of the Onsager reciprocal relation in transport constitutive laws of the Stefan-Maxwell form, 13 emphasizing that no particular choice of excess fluxes needs be made when developing constitutive laws for transport. In the discussion of momentum it will be convenient to introduce a particular set of excess fluxes. We will use the set of molar fluxes relative to the mass-average velocity, which are constrained by a kinematic relation but wish to emphasize that this choice is purely one of convenience. In particular, these excess fluxes preserve the standard definition of the substantial derivative, allowing one to write the familiar expression for the balance of species i in a frame following the bulk motion of the material.

Momentum Balances
An overall local momentum balance accounts for how external forces drive material motion. When deriving this balance, which is typically expressed for pure phases by the Cauchy equation, it is not so clear what convective velocity should be involved. (This issue received substantial attention from Brenner in the early 2000s. 27 ) By starting with species material balances above, the mass-average velocity v was successfully identified as the relevant velocity in the traditional masscontinuity equation. A similar tack will be taken here, and the total momentum balance will be derived by combining individual species momentum balances, which avoids giving precedence to a particular bulk velocity.
Thermodynamic force.-One must first consider how thermodynamic gradients manifest as forces, i.e., how such gradients locally generate momentum. An apparent body force arises from spatial gradients in the available mechanical energy G mech , and species-specific thermodynamic driving forces are produced from local gradients in the rest of the available energy (G − G mech ).
One can divide Euler's equation for the grand canonical ensemble (Equation 20) by m T , take a total differential, and insert Gibbs-Duhem Equation 45 to get an expression for ρdĜ. Considering only the parts that vary when temperature and composition are held constant, one identifies ρ ∇Ĝ mech = ∇ p + : ∇ τ. [83] Dividing both sides of this expression by ρ yields a formula for local body acceleration, and multiplication by ρ i delivers the mechanical force per unit volume that acts on species i. Since ρ i /ρ is the mass fraction ω i , the mechanical force per unit of volume that acts on species i must be where the latter two terms in the brackets arise from a tensor identity. Thermodynamic forces are identified via the procedure of Hirschfelder, Curtiss, and Bird. 29 The spatial variations of intensive properties in each direction are taken to satisfy the differential relationship established by Gibbs-Duhem Equation 45. After introducing the Euler equation for entropy from Equation 32, then dividing both sides by V , one can identify a thermodynamic force d i with units of force per volume, [85] This definition recasts the Gibbs-Duhem equation as The local balance of thermodynamic forces on individual species can be viewed as a statement of Newton's third law of motion for systems in microscopic equilibrium. 30 Momentum continuity.-As well as being the mass flux of i, ρ i v i = m i N i represents the contribution that species i makes to the local momentum density. b Momentum continuity for an individual species in a phase can then be expressed as Brenner's work suggests that this statement amounts to a physical postulate, which may require additional justification. 27 Here we have left the microscopic definition of v i ambiguous intentionally, with the understanding that an appropriate choice might address Brenner's concerns. For example, choosing a microscopic averaging process such that v i is, in Brenner's terms, a 'volume velocity' of species i may resolve the problem. Of course, this contradicts the usual interpretation of v i as a number-average velocity over all particles of type i.

E3656
Journal of The Electrochemical Society, 164 (11) E3647-E3660 (2017) Cauchy stress is apportioned among species here in the same way as it was when defining the diffusion driving force in Equation 85. The vector b i represents all the body forces imposed on i by auxiliary fields not already considered in the Thermodynamics section. Examples of such body forces include gravity and the Lorentz force, both of which are governed by complementary sets of field equations that can be appended to the present model without contradiction. Summation of Equation 87 over all species yields an overall momentum balance. Unlike the summation of species material balances, which results in the familiar continuity equation, the overall balance derived from species momenta contains an additional convective term, because an identity that can be derived by applying Equations 80 and 81. The diffusion stress , defined as accounts for inertia due to excess momentum convection, which exists because species i carries momentum with itself at v i , rather than at the bulk velocity v. c In the overall momentum balance, Equation 86 shows that the contributions of the thermodynamic forces vanish; application of mass continuity from Equation 71, introduction of the convective derivative from Equation 72, and insertion of Equation 2 further simplify the sum to This looks much like the familiar Cauchy equation, save for the added diffusion stress.

Energy Balances
The heat equation ultimately arises as a consequence of the balance of total energy, which is the sum of specific kinetic energy and specific internal energy. Total specific energyÊ can accumulate locally within a material element as a consequence of heat exchange with its surroundings, exchange of mechanical energy with its surroundings, or energy absorption due to the action of external fields: This equation introduces the heat flux q, whose definition will be refined somewhat in the sequel. Equation 91 basically matches the forms used by Hirscshfelder et al., Bird et al.,and Gyarmati,29,31,32 except that it includes energy gained as a consequence of diffusion stress. Continuity of kinetic energy is a consequence of the modified Cauchy equation. Dot v into both sides of Equation 90, then identify the specific kinetic energyK to show that Again this expression matches canonical results, 31  yields a balance of internal energy, which has been simplified with the tensor identity The form of the last term on the right of Equation 95 can be clarified in several ways. To emphasize that the antisymmetric part of ∇ v does not contribute to the double dot product, replace ∇ v with the rate-of-deformation tensor a quantity familiar from fluid mechanics. In place of the Cauchy stress, introduce pressure and deformation stress with Equation 2, then apply the identity I : ∇ v = ∇ · v to simplify the pressure term. Finally, put the diffusion stress back in terms of fluxes using Equation 89 and apply tensor identities to obtain It only remains to express the left side of Equation 97 in terms of thermodynamic basis variables and quantities that describe the elastic state. Following Hirschfelder et al., use Equation 59 to express a relationship among the substantial derivatives of the quantities involved, i.e., [98] The derivatives of ρ and c T y i can be eliminated in favor of flux divergences using Equations 73 and 82, respectively. With this in hand, one can eliminate DÛ /Dt with Equation 97, after which some straightforward rearrangement gives [99] A final simplification is provided by recognizing that the heat flux q has contributions both from irreversible heat flux -by mechanisms like conduction -and the latent heat carried by excess species fluxes.
To separate these, introduce the irreversible heat flux which isolates the flow of heat associated with purely dynamical processes. (Hirschfelder et al. introduced q ; 29 whereas q is the excess heat flux in excess of the enthalpy convected at v, q is the excess flux relative to convection at v h .) Putting Equations 99 and 100 together, one achieves a thermal energy balance applicable to elastic, multicomponent materials: [101] This form of the general heat equation matches the expression given by Bird, Stewart, and Lightfoot, 31 [102] When dotted with J i and summed over all species, the term at the far right vanishes because of kinematic relation 81. Thus an equivalent form of the heat equation is [103] Reading in order from left to right, the terms signify: heat accumulation, heat efflux due to dynamical processes, reversible energy accumulation associated with compression and deformation (that is, entropic heat from mechanical processes), energy dissipation from diffusional stress, diffusion, and friction, and heating from excess entropy convection.

The Dissipation Function
The development now culminates with a transient statement of the second law of thermodynamics. Within Onsager's transport theory, 19 and within the field of nonequilibrium thermodynamics in general, 34 this law is derived by performing a dynamical entropy balance, then using thermodynamic relationships and the continuity equations for mass, momentum, and energy to express the entropy generation in terms of a sum over pairs of conjugate fluxes and driving forces.
A general continuity equation for entropy can be written as Another formula for the substantial derivative of entropy can be obtained by combining the substantial derivative of Equation 57 with Equations 71 and 82 to express DÛ /Dt, which can be eliminated with a combination of Equations 97 and 100 to get [106] Elimination of DŜ/Dt between Equations 105 and 106 and insertion of Equation 102 yields [107] One final simplification allows the energy dissipation to be expressed in terms of independent fluxes and driving forces: application of kinematic relation 81 and Gibbs-Duhem Equation 86 can be used to eliminate the n th term from the sum in Equation 107. Ultimately the energy dissipation is given as a sum over pairs of independent fluxes and their conjugate driving forces: In this expression the first term on the right illustrates the fundamental fluxes and forces that appear in the constitutive laws for Onsager-Stefan-Maxwell diffusion; 15 the term with q underpins Fourier's law of heat conduction, and can be used in combination with the diffusion term to rationalize the Soret/Dufour effects; 2,12,29,31 the last term justifies the Navier-Stokes constitutive law for Newtonian fluids. Several new quantities involving the deformation strain appear in the equation, providing a logical basis for the formulation of dynamical constitutive laws that apply to viscoelastic or plastic materials. The term with ˙ · m i J i is a cause for concern because it suggests that diffusion stress, as well as the thermodynamic force d i , can drive mass diffusion. This factor is typically not included in the Onsager-Stefan-Maxwell formalism. Its addition may lead to difficulties, since it will result in transport constitutive laws that mix fluxes and driving forces, which Monroe and Newman have shown can affect the symmetric Onsager reciprocal relation one expects among the transport coefficients. 13 Probably this is not an issue, however, because the derivation of a reciprocal relation via Onsager-Casimir fluctuation theory involves linearizing the transport constitutive laws around an equilibrium state. 15 For equilibrium fluctuations ˙ · m i J i is a term of higher than linear order.
Another way of handling the problem with the ˙ · m i J i term is to introduce the diffusion stress, and rewrite Equation 108 as [109] Here diffusion stress is identified as a type of momentum flux, with the rate-of-deformation tensor providing the conjugate driving force. This formulation seems to require a separate constitutive relationship to connect with ˙ , however.
In closing, the second law of thermodynamics can be stated simply as with equality holding only at equilibrium. This completes the development of the general equilibrium governing equations and dynamical balance laws applicable to elastically deformable multicomponent electrolytes.

Conclusions
Most of the quantities introduced here are familiar from prior implementations of multicomponent diffusion theory. [29][30][31][32][33][34] The new model essentially reduces to those prior developments with only a few minor modifications, most of which involve the strain.
The development here has provided the groundwork needed to state a sufficient set of balances and variable definitions for a description of simultaneous transport of material, momentum, and heat within an electrochemically active material. For an n-ary material in a ddimensional geometry, Equations 15, 28, 55, 73, 81, 89, 90, 96, and 103, along with n − 1 of Equations 82 and n of Equation 85, constitute a set of 4 + (n + 2d)(1 + d) scalar equations (each vector equation being counted d times, and each tensor equation, d 2 times); these involve 4 + (n + 2d)(1 + 2d) scalar dependent variables (with vector quantities counted d times, and tensors, d 2 times): {y i } n , T , p, c T , ρ, τ, v, { J i } n , q , { d i } n , , ˙ , and . The electrical state of the material is described by adopting Equations 64 and 65, which introduce the excess charge density ρ e and current density i without changing the balance of equations and unknowns. In some cases it may be helpful to describe kinematics with the volume-average velocity, which can be computed by appending Equation 62; this allows one species material balance to be replaced with Equation 63.
The governing equations also involve various equilibrium material properties, which, being based in the Gibbs function, depend at most on temperature, stress, and composition: Description of the electrical state also involves the species charges {z i } n . If a volume balance is adopted, one may also require the material's modulus K . One should be careful to ensure that the electrochemical potentials relate to partial molar volumes and entropies through appropriate Maxwell relations (Equations 24 and 31). The Euler equation for V , Equation 25, shows that α V and K , as defined by Equation 30, are determined by the constitutive laws for electrochemical potential; the Euler equation for S, Equation 32, leads to a similar connection between the set of μ i andĈ p . All the partial molar properties must be properly constrained by Gibbs-Duhem relations (Equations 45, 46, and 49); these will in turn ensure that the diffusion driving forces satisfy Equation 86, which is a fundamental law in Onsager-Stefan-Maxwell transport theory. 19 Finally it is possible to provide a count of how many flux constitutive laws are needed to close the governing system of transport equations. The number of scalar dependent-variable components is greater than the number of scalar governing equations by nd + 2d 2 . Of these remaining degrees of freedom, d 2 can be specified by adopting an equilibrium constitutive law for deformation strain like Equation 41. Inspection of the entropy generation in the form of Equation 108 suggests that transport constitutive laws should close nd degrees of freedom by relating fluxes q and { J i /c T y i − J n /c T y n } n−1 to the independent forces − ∇ ln T and { d i + ˙ · m i J i } n−1 ; the last d 2 degrees of freedom are closed by expressing the flux τ as a function of the rate-of-deformation tensor ˙ , the substantial derivative of the deformation strain tensor , and ∇ · v. Explicit formulation of these thermodynamic and transport constitutive laws will be a topic of future work. of the University of Michigan. The authors thank Prof. C. W. MacMinn (Oxford, Department of Engineering Science) for his helpful early insights about mechanical principles, and Prof. C. P. Please (Oxford, Department of Mathematics) for his suggestion to investigate species momentum balances, which led us to identify the diffusion stress.

Appendix: Mathematical Properties of Tensors
The 3-tuple { e i } 3 represents a right-handed, orthornormal set of basis vectors. Between a vector v and a tensor σ there are two ways of forming a new vector with a scalar product, [A1] The transpose of σ is the unique tensor σ T that satisfies for every v. Between tensors σ and τ, forms a tensor; this product does not generally commute.
The trace of tensor σ is defined as Generally tr( σ) = tr( σ T ) and tr( σ · τ) = tr( τ · σ). The adjugate trace of σ is defined in terms of the trace and scalar product as and the determinant of σ is The determinant distributes over the scalar product of tensors as det( σ · τ) = det( σ)det( τ).
The identity tensor I satisfies v · I = I · v = v for all v and σ · I = I · σ = σ for all σ. In component form, where δ i j is the Kronecker delta, equal to 1 if i = j and 0 otherwise. The inverse of σ is the unique tensor σ −1 for which Note that tr[adj( I)] = 3 and det( I) = 1. Since det( σ · σ −1 ) = det( σ)det( σ −1 ) = 1, the determinant of a tensor must be nonzero for its inverse to exist.
The vertical scalar product between tensors σ and τ is given by With this definition one can also write the trace of σ as a vertical scalar product with the identity, tr( σ) = σ : I , whence tr( I ) = I : I = 3. Properties of the trace imply that the vertical scalar product commutes, σ : τ = τ : σ.
One can understand how the components of a tensor transform under a rigid rotation of coordinate axes by observing that the act of rotating a vector v without changing its magnitude can be expressed by projecting a rotation tensor Q onto v, where Q T = Q −1 and det( Q) = 1. The most primitive rotation is expressed by the identity tensor I , which merely maps a vector or tensor into itself.
If T is any tensor, then rotation of the basis vectors under T by Q sends T into a new tensor T according to [A10] whose components stand in the rotated coordinate frame. The trace, adjugate trace, and determinant defined by Equations A4 through A6 remain unchanged under transformations like the one in Equation A10, and are therefore known as a tensor's three principal invariants, which express its essential physical content. The spectral theorem from linear algebra codifies a few general consequences of tensor symmetry: given any tensor T = T T , with components over a right-handed, ordered orthonormal basis { e 1 , e 2 , e 3 }, there is some rotation tensor Q T with det( Q T ) = 1 such that The ordered list { e T 1 , e T 2 , e T 3 }, where e T i = Q T · e i , represents a rigidly rotated, righthanded, ordered orthonormal set of principal axes for tensor T , and the nontrivial entries The principal stresses of a tensor T generally relate to its three principal invariants through [A12]