Incremental least action principle in the framework of thermodynamics of irreversible processes

In this paper, an incremental least action principle is proposed using the framework of the thermodynamics of irreversible processes. First, we establish that an extensive thermodynamic potential deﬁned over an inﬁnitesimal volume satisﬁes a differential conservation law with close similarities to the Liouville theorem but extended to irreversible processes. This property is applied to a generalized thermodynamic potential depending on equilibrium variables and nonequilibrium ﬂows. It allows the formulation of an absolute integral invariant (AII) that is shown to have a broader ﬁeld of application than the Poincaré-Cartan integral invariant of dynamic systems. Once integrated over a ﬁnite volume, it naturally deﬁnes an integral functional that fulﬁlls an incremental least action principle. The Fréchet derivative of the Euler-Lagrange equations associated with the functional is calculated, and its self-adjointness is shown to be equivalent to the symmetry of the classical Tisza and Onsager matrices which link respectively extensive variables to intensive variables and nonequilibrium ﬂows to generalized forces. Finally, the proposed AII and least action principle are formulated for the case of a simple physical process (heat conduction), to illustrate (i) its link with the extended irreversible thermodynamics and (ii) its applications to numerical simulations.


I. INTRODUCTION
Despite the warning of the French mathematician Henri Poincaré, who stated in 1908 that "irreversible phenomena cannot be explained by means of Lagrange equations" [1], it has been the dream of many physicists to reconcile the thermodynamics of irreversible processes with the least action principle.Beyond its philosophical foundation, the use of an integral functional-whose extremals coincide with the governing equations of a physical system-turns out to have important practical interests like the method developed by the Swiss mathematician Walther Ritz and known as the "Ritz method" [2], valid when the second variation of the functional is positive definite.In the mechanics of elastic bodies, this method consists of minimizing the elastic energy for which a simplified guess of the displacement field has been introduced.It leads to the identification of an approximate displacement field with good precision and avoids the calculation of the true solution to the equilibrium equation (see, e.g., [3]).In the mechanics of punctual masses, calculating the Lagrangian (as the difference between kinetic energy and potential energy), and writing the associated Euler-Lagrange equations have proven to be very useful in efficiently obtaining the equations governing the motion of these masses [4].Unfortunately, not all systems of partial differential equations (P) correspond to a system of Euler-Lagrange equations associated with a Lagrangian L. This possibility actually relies on the self-adjointness of the system (P).Classically, we observe that the differential operator appearing in the equations governing the evolution of a system having an irreversible evolution is not self-adjoint.More precisely, the theorem presented in [5] and given in a more general form in [6] states that a Lagrangian L exists for a system (P) if and only if its Fréchet derivative D P is self-adjoint, i.e., equal to its adjoint D * P .We recall that the adjoint D * P of a linear differential operator D P having constant coefficients is obtained by switching the sign of odd-order derivative operators and by transposing all the coefficients appearing before derivative operators.For instance, when equations contain first-order derivatives in time, irreversibility of time breaks the symmetry t → −t and the Fréchet derivative D P cannot be self-adjoint.As an illustrative example, we consider here the simple case of an oscillator of mass m with only one degree of freedom q(t ), linked to a spring of stiffness k and put inside a fluid bath involving a viscous friction λ.The differential equations satisfied by q(t ) and its Fréchet derivatives are given in the first, second, and third rows of Table I, respectively ("irreversible case").It can be seen that a Lagrangian cannot exist for this problem since D P = D * P , unless λ = 0, which corresponds to the reversible case with no friction.Many approaches have been proposed in the literature to formulate TABLE I. Comparison of the three main equivalent Lagrangian approaches found in the literature, for a linear damped (λ = 0) oscillator having only one degree of freedom.For comparison, the nondissipative reversible case (λ = 0) is introduced.The descriptions of the three approaches are as follows.(a) Integrating factor (IF) strategy: Create an explicit time dependence.(b) Adjoint system (AS) strategy: Introduce a new variable u(t ) such that u(τ − t ) = q(t ).(c) Rayleigh function (RF) strategy: Reversible Lagrangian supplemented by a Rayleigh function R( q) (not a true variational principle).
Reversible Case (λ = 0) Irreversible Case (λ = 0) Comment Governing equation P = m q + kq = 0 P = m q + λ q + kq = 0 D tt , D t : Total time derivative operators Fréchet derivative Rayleigh function (RF) strategy q(t ) E q (L RF ) = ∂R ∂ q "equivalent Lagrangian formulations" of this problem using the property that a Lagrangian L associated with equations (P ) is said to be equivalent to a Lagrangian L associated with (P) if the solutions to the initial problem (P) are equivalent to those of (P ) or included in them.In the literature, one can distinguish three main strategies for obtaining an equivalent Lagrangian formulation.The first one consists of multiplying the initial system (P) by a nondegenerate matrix such that the new system is self-adjoint.An example of this strategy has been proposed in [7] for the oscillator previously presented [see Table I, irreversible case (a)], but also in thermics for the Cattanéo-Vernotte equation and for an electrical circuit containing resistors [8].In the latter case, the integrating factor introduced for the multiplication can be expressed as a matrix exponential.The main disadvantage of this strategy is that it is limited to regular systems of partial differential equations [5], and it cannot be applied to diffusion processes which are associated with a differential operator having only a first-order derivative in time.The second strategy consists of introducing new variables called "dual variables" and satisfying an adjoint system (Q) such that the entire system (P, Q) is self-adjoint.Methods can be found in the literature to calculate the adjoint system [9,10].This approach has been used in thermodynamics through the work of Anthony [11] who introduced a thermion field χ and a phase field ξ such that φ = √ T e iξ , where T is temperature.Concerning our oscillator example, a "bi-Lagrangian" has been proposed in [12] and is presented in Table I, as the irreversible case (b).The main concern that arises from this strategy is the physical meaning of the dual variables, which have an evolution inverse to time.In the case of the present oscillator example, introducing t = τ − t in u(τ − t ) = q(t ) leads to u(t ) = q(τ − t ) and clearly highlights this reverse time evolution.Since the oscillator is damped by the fluid, q(t ) exponentially decays toward zero and the quantity u(t ) exponentially increases with time, which has no trivial physical meaning.The last strategy that can be found in the literature relies on dissipation potentials (sometimes called Rayleigh functions).It consists of identifying the first variation of the action integral S = L dt with the increment of energy dissipated by the system integrated over dt.It thus introduces the missing irreversible term in the Euler-Lagrange equations.To illustrate this approach, we have presented the Rayleigh function of our oscillator in Table I, as irreversible case (c).This strategy has also been used in continuum mechanics [13][14][15][16], and has an incremental formulation in terms of differential geometry [17].However, we point out that despite its great efficiency in practice, this formulation is not a true variational principle in the sense that the Euler-Lagrange equation supplemented by the irreversible term derived from the Rayleigh function does not correspond (a priori) to the first variation of any integral functional.
If we look at these three equivalent Lagrangian formulations, we observe that the self-adjointness condition has always been circumvented by mathematical tricks with nonphysical characteristics or away from the initial idea of having a variational problem.In our opinion, this is mainly due to the fact that the physical meaning of self-adjointness remains not fully understood yet.One of the aims of this paper is to shed some light on this property.To achieve this objective, four main steps are identified.At first, we establish that every thermodynamic potential defined over an infinitesimal volume satisfies a differential conservation law as soon as it fulfills the axioms of Callen's approach [18].It is shown that this step requires a new definition of extensivity, which is valid when extensive variables become differential forms.Second, we apply this differential conservation law to a generalized potential in the sense of extended irreversible thermodynamics, and formulate a general absolute integral invariant (AII) in the framework of the thermodynamics of irreversible processes.This AII is shown to have a broader field of application than the Poincaré-Cartan AII (formulated for reversible dynamics) since it is valid for irreversible processes.In a third step, it is demonstrated that our AII can be integrated over a finite volume to obtain an incremental least action principle formulated as an integral functional T .The Fréchet derivative of the Euler-Lagrange equations associated with the stationarity of T is calculated to propose a physical meaning of self-adjointness.Finally, a detailed example corresponding to a simple physical process (heat conduction) is considered to illustrate the approach, clarify technical aspects, and discuss the link of the present approach with the extended irreversible thermodynamics [19].

A. A differential form of extensivity
We consider here a thermodynamic system occupying a material elementary volume associated with the differential 3-form = dx ∧ dy ∧ dz and remaining at equilibrium.In the definition of , d denotes the exterior differential and ∧ the exterior product of the differential forms.In accordance with Callen's axiomatic approach of thermodynamics [18], can be described by a set of n independent extensive variables, and a thermodynamic potential S(E ) (entropy or energy) containing all the information on the system and satisfying the extensivity property, i.e., the proportionality to the amount of matter involved: If we introduce the specific potential s and the specific coun- and (ii) if the following property is satisfied: in which L[X α ]S stands for the Lie derivative of the differential form S with respect to the vector field X α .Throughout this paper, if not otherwise stated and as in Eq. ( 4), the Einstein summation convention on repeated indexes is adopted.In Eq. ( 4), I k are the dual intensive variables of E k in a thermodynamic sense and will be defined more precisely in the next section.Equation ( 5) is simply a differential formulation of the Euler extensivity property: equivalent to (2) with λ = e μα , since X α is the generator of the Lie group corresponding to multiplication by some factor e μα , μ being the parameter of the group.The writing of ( 5) only assumes that the extensivity is defined by multiplying e k by e μα at constant volume , which is equivalent to multiplying E k by the same factor.

B. Differential conservation law
The aim of this section is to show that the previous redefinition of extensivity [Eq.( 5)] preserves the main results of thermodynamics defined for finite volumes.However it also shows an interesting local invariance property.First we demonstrate that S satisfies the fundamental Euler equation: The proof of this result is given in Appendix A. It shows that the classical expression of I k as derivatives of s with respect to e k [18] is preserved.Second, from Eqs. ( 4) and ( 7), one obtains Subtracting dS from both sides of the equation completes the proof of the Gibbs-Duhem relation in differential form: Note that introducing E k = e k in ( 9) and collecting terms in leads to the classical Gibbs-Duhem equation e k dI k = 0 as formulated in [18].The present Eq. ( 9) is a more general formulation since it accounts for the fact that E k is defined over the infinitesimal volume .Third, taking the exterior differential of (4) leads to We qualify this last relation as the "Maxwell-Liouville" relation.Indeed, its physical meaning reveals that as long as S remains a thermodynamic potential (with extensive arguments) and satisfies the Maxwell relations (d 2 S = 0; see Appendix B), then the 5-form dI k ∧ dE k (interpreted as the elementary volume of the space {E k , I k }) is preserved and remains equal to zero.We emphasize that Eq. ( 10) is valid for any infinitesimal volume , and will be further integrated over a finite volume to define a functional satisfying an incremental least action principle.In that sense, the quantity dI k ∧ de k appearing before in (10) can be viewed as a specific Lagrangian as introduced in classical variational calculus [6].One notes the high similarity with the Liouville theorem which states that the volume of the phase space d p ∧ dq is preserved for any frictionless movement associated with the position q and the momentum p of a punctual mass [17].Finally, we see that the extensivity of S seems to be unnecessary to obtain Eq. ( 10), but is, in fact, essential to the formulation of Eq. ( 10) for any Legendre transform of S (see Sec. III C).For the sake of clarity, Fig. 1 gives a geometrical interpretation of Eq. ( 10) when the potential S depends only on two extensive variables E 1 and E 2 .The evolution of the potential 1. Geometrical interpretation of Eq. ( 10): when the potential S depends only on two extensive variables: Note that due to the antisymmetry of the exterior product, surfaces may be understood in an algebraic sense.Consequently, in the present graphs I 1 vs E 1 and I 2 vs E 2 , blue angled hatching surfaces (resp.pink filled surfaces) must be equal, following the Gibbs-Duhem Eq. ( 9) [resp.Eq. ( 10)].equation, Eq. ( 7), is represented as the contributions to the potential S from the red areas: S 1 and S 2 .The Gibbs-Duhem equation [Eq.( 9)] is shown as the property that the sum of the blue areas has to be equal to zero when considering surfaces in the algebraic sense: 10) is illustrated as the balance of the purple areas:

A. Extended thermodynamic potentials
The results of the previous section are now applied to the thermodynamics of irreversible processes.The term "extended thermodynamic potentials" in the title of the section is borrowed from the "extended irreversible thermodynamic" approach of Jou et al. [19].As previously introduced, we consider the thermodynamic system occupying the infinitesimal volume .We temporarily assume that the system is at equilibrium, or near equilibrium.We have seen previously that this system can be described by the variables {E k , k ∈ [[1; n]]} and the potential S(E ).Despite the generality of the Callen axiomatic approach, the thermodynamic potential S(E ) cannot describe dissipative changes of .Indeed, S(E ) does not contain any information about timescales for which the system would return back to equilibrium if it were driven away from it.It is thus necessary to introduce m supplementary variables ]} to describe the irreversible evolution of [20].In their foundation of extended irreversible thermodynamics (EIT), Jou et al. [19] generalize the notion of entropy outside of equilibrium.To achieve this, they introduce an entropy term S EIT that depends on a new set of variables (E, J) (thus extending the initial set E) as follows: In the present study, this idea is retained but in a slightly different manner.Indeed, in the same spirit as Jou et al., the set of extensive variables used as arguments of the generalized potential S is split into two sets of variables (E, ): However, here, the new variables = { α , α ∈ [[1; m]]} are defined as cumulated flows with respect to time (see their specific definition below).These kinds of variables have actually already been used in [21] to transpose a thermoelastic problem in variational form.Interestingly, some Tonti charts [22] exhibit unnamed variables which could coincide with our cumulated flows .This is the case for the integral of the temperature gradient with respect to time in the case of thermal conduction (see Table 13.4 in [22]).We introduce the specific cumulated flows ψ (defined by α = ψ α ) and the specific flows (defined by J α = j α ), linked through the differential equation Using the set of variables of Eq. ( 12), and introducing the specific generalized potential s defined by S = s(e, ψ) , it is possible to develop the second exterior differential of S as with Further development of Eq. ( 14) leads to since dψ α = j α dt.In accordance with the Maxwell-Liouville equation [Eq.(10)] reconsidered here with the extended set of variables (I k , e k ) ≡ (I k , e k , X α , ψ α ), the quantity vanishes as long as S satisfies the Maxwell relations d 2 S = 0.A particular formulation of this property involving the Tisza and Onsager matrices is presented in the next section.

B. An absolute integral invariant
In this section, we analyze the particular but frequent case for which the specific potential s is split into an equilibrium potential s eq that depends only on e, and an irreversible term s irr that depends only on ψ: s(e, ψ) = s eq (e) + s irr (ψ). ( The time evolution of the dual variables I is given by the differentiation of the first equation of ( 15) accounting for Eq. ( 17) and multiplied by : which, using (B2), is strictly equivalent to In Eq. ( 19), the generalized stiffness matrix T = [T k j ] as defined by Tisza is introduced [23].This matrix is necessarily (i) symmetric while the equilibrium potential s eq (e) satisfies the Maxwell relations and (ii) positive definite as long as the equilibrium around which the system is considered is stable.It is commonly accepted in the thermodynamics of irreversible processes (see, e.g., [20]) that the nonequilibrium forces X α and the nonequilibrium flows j α are linked via a coupling matrix L = [L αβ ].As for T , L has to be symmetric to fulfill Onsager reciprocity relations and positive definite to ensure the validity of the second law of thermodynamics.We consider here the incremental writing Based on Eqs.(19) and (20), it is possible to describe the thermodynamic evolution of the system by introducing a vector field ζ defined by Indeed, the exponential of this vector field is given by the integration of the differential system [24]: which contains Eqs. ( 19) and ( 20) but accounts for any evolution of e and X .The main benefit of introducing ζ is that the 5-form defined by Eq. ( 16) and rewritten here for convenience, is an AII for the vector field (21), i.e., L[ζ ] = 0, if and only if T and L are symmetric.The proof of this theorem is given in Appendix C. It gives a more straightforward physical meaning of the Maxwell-Liouville relation (10) in the case of Eq. ( 17) since the symmetry of T and L are directly linked to the Callen axiomatic approach and the Onsager relations.

C. Invariance by Legendre transform
One could object that the previous result is only valid for thermodynamic potentials that depend solely on extensive variables.In fact, the formulation of is independent of the state variables and the nature of the chosen thermodynamic potential (Legendre transforms of energy, or Legendre transforms of entropy, also called Massieu functions).Indeed, can be viewed as the second exterior differential of the generalized potential S = s(e, ψ) but also as the second exterior differential of its full Legendre transform F = f (I, X ) .The relations defining e and ψ are then given by from which one can obtain This result remains actually true for any partial Legendre transform of S. In practice, this invariance property widely extends the field of application of the Maxwell-Liouville theorem.In particular, the integral of the differential 5-form over the volume has important implications for the formulation of convex integral functionals since a thermodynamic potential is convex with respect to its extensive parameters and concave with respect to its intensive parameters (see Sec. III D).

D. An incremental least action principle
There is a very close link between the Maxwell-Liouville theorem and variational calculus.To briefly discuss this point and for the sake of clarity, it is assumed in this section (but without loss of generality) that the specific flows j and the generalized thermodynamic forces X are vectorial quantities and that the number of extensive variables E k equals the number of flows, i.e., n = m.The index notations will thus be changed since the same indexes can be used to loop over the variables (I, e) and (X , j).We also assume that the system is now occupying a finite volume associated with the infinitesimal volume dω, and subsequently that each thermodynamic variable becomes a field depending on the space variables x, y, z.Starting from Eq. ( 23), it is possible to associate an integral functional T with , defined by Another form of T can also be introduced using flows calculated at the end of the time increment δt as Both functionals may be supplemented by boundary terms ∂ • • • da not written in Eqs. ( 26) and ( 27), for readability.For many dissipative processes, the thermodynamic forces are defined as the gradient of intensive variables.Accordingly, we assumed here that indexes i ∈ [[1; n]] are chosen such that δX i can be rewritten as If T is now considered as a functional of δI only, one obtains If δI is split into two parts, with δI ∞ satisfying the Euler-Lagrange equations associated with the Lagrangian of T , it can be shown by integration by parts that in which W (δI * ) is the new functional defined by If the equilibrium around which the system evolves is stable, the Tisza matrix T and its inverse are positive definite.The coupling matrix L is also positive definite to account for the second law of thermodynamics, and we have W (δI * ) 0, ∀δI * .If the increments δI are chosen in a set U such that the boundary terms of (32) vanish, it can be concluded that and that T is minimal when δI = δI ∞ .Interestingly, Eq. ( 31) is simply an implicit (in a numerical sense) writing of the balance equation associated with the specific quantity e i during δt and can be rewritten . Consequently, the variational principle defined by Eq. ( 29) can be viewed as an incremental least action principle.It has a practical interest if each field δI i , as for the Ritz method, is developed as a linear combination of some functions N ik (x, y, z): γ k being unknown coefficients.Inserting (36) into (29) and minimizing the result with respect to γ k lead to a linear system in γ k .Putting back the values of γ k in (36) gives an approximation of δI ∞ i .This strategy, repeated for each time increment δt and by updating I with δI, will be qualified as the "incremental Ritz method" (IRM) hereafter.The complete study of the functionals ( 26) and ( 27) and its possible link with the extremal problems developed in the work of Sewell [25] goes beyond the scope of this paper.However, a detailed formulation of T and its subsequent minimum incremental principle will be proposed for the heat equation in Secs.IV B and IV C, as a detailed example.

E. The physical meaning of self-adjointness
A direct calculation of the Fréchet derivative of Eq. ( 31) using the formula given in [6] leads to in which (i) D k denotes the total derivative operator with respect to the space variable x k (using the new notation x 1 = x, x 2 = y, x 3 = z), and (ii) D J denotes the total derivative operator with respect to a multi-index J (e.g., ).We also recall that for any vector component The adjoint of D can be calculated applying the direct formula of the same author [6]: with (−D) J = (−1) s D J , s being the number of indexes in J.
The comparison between Eqs. (37) and (38) shows that the Fréchet derivative of = { i , i ∈ [[1; n]]} [given by Eq. ( 31)] is self-adjoint if and only if which corresponds to the symmetry of the Tisza and coupling Onsager matrices.This result may appear trivial since i corresponds to the Euler-Lagrange equation of the functional T with respect to the intensive variable δI i .However, it lends an interesting physical meaning to the self-adjointness property.

F. Link with the Poincaré-Cartan invariant
The 5-form may be rewritten in terms of dissipation if a linear relation between flows and forces is considered, as in [19].Indeed, if Eq. ( 23) is rewritten accounting for j α = L αβ X β (with L a constant, symmetric, and positive-definite matrix), one obtains Introducing the specific dissipation as follows, Comparison between the absolute integral invariant formulated as if it were a Massieu function and its EIT counterpart in the case of thermal diffusion [19].In this table, u denotes the internal specific energy, v the specific volume, q the conductive thermal flow, T the temperature, κ the thermal conductivity, τ 1 some time constant, S the entropy, and t the time.The relation q = −κ∇T has been used in the formulation of .

Framework
Potential or Invariant Equilibrium Term Dissipative Term EIT, Eq. (2.40) of [19] since d = 0.This last form exhibits a high similarity with the Poincaré-Cartan AII [26]: One of the main differences between these two AIIs lies in the fact that is homogeneous to an energy measured in J (since is in J m −3 s −1 ) while PC is homogeneous to an action measured in J s.However, a close similarity can be emphasized in the case of mechanics for a punctual mass submitted to potential forces only (no dissipation, i.e., dH = 0 and d = 0).In that case, the only extensive quantity E 1 = e 1 is the momentum p in kg m s −1 and its thermodynamic intensive dual is the velocity v in m s −1 .The analogy between PC and is underlined by the equivalence Since v = q in this example, and since dH ≡ d ( ) by comparing (43) and ( 44), it all happens as if would be the time derivative of PC .One notices that the Poincaré-Cartan AII was initially written in the case of potential force fields, i.e., without any dissipation.The major contribution of the present approach is that the integral invariant does not require any assumption of reversibility.In that sense, its field of application can be viewed as more general than that of PC .

A. Similarity to the extended irreversible thermodynamics
To discuss the similarity between the present approach and extended irreversible thermodynamics [19], the absolute integral invariant is here written in the case of heat diffusion as if it were a Massieu function.Below and in Table II, e 1 = u denotes the internal specific energy, v the specific volume, j 1 = q the conductive thermal flow, T the temperature, κ the thermal conductivity, X 1 = ∇T −1 is the thermodynamic force associated with q, and the thermodynamic dual of the specific energy u is the quantity I 1 = 1/T .The relation q = −κ∇T has been used in the formulation of (see Table II) while it can have a more general form as discussed in [19].Consequently, accounting for ∇T −1 = −∇T /T 2 = q/κT 2 , one can rewrite Eq. ( 23) as Interestingly, a similarity between and its EIT counterpart d 2 S EIT (see Table II) can be highlighted.First, the equilibrium term of is simply given by eq = d and coincides with the second exterior differential of entropy S eq as it is classically defined at equilibrium As concerns the dissipative part of , shown in the last column of Table II, other similarities emerge: τ 1 ≡ dt, v ≡ , and (q • q)/2 ≡ q • dq.The main difference between these two dissipative terms is that in the present approach, no time constant τ 1 needs to be introduced.This intuitive comparison between and d 2 S EIT in the case of thermal conduction could be more deeply analyzed and generalized to other physical processes such as mass flow, or poromechanics.

B. An incremental minimum principle
In this section, we turn back to an energy description by considering I 1 = T , e 1 = s (specific entropy), j 1 = q, X 1 = g/T with g = −∇T , and the specific heat c p such that One considers that the system is submitted to a given flow q d over a part ∂ q of its boundary ∂ and a given temperature T d over ∂ T with the conditions One introduces a time discretization and the notation δ • will be used to denote the increment of any quantity • over the time increment δt.In the same manner, (•) k refers to any quantity evaluated at time t k .We rewrite Eq. ( 46) in an energy form rather than in an entropy form ( ≡ d 2 S while ≡ d 2 U with U the internal energy generalized as in [19]), i.e., Following from Eq. ( 29), it is possible to define the integral functional for known T k and q k (temperature field and thermal flow at time t k ): and the set of admissible temperature increment fields δT k satisfying the boundary conditions over the time increment δt: where the superscript "ad" stands for "admissible".As in Sec.III D, it can be shown by integration by parts that for any element of U ad k , T can be rewritten as which is the entropy balance during the time increment δt.
Since T k > 0 over , and since κ k > 0 (second law of thermodynamics) and c p > 0 (stability of equilibrium), it can be concluded from Eq. (56) that T (δT k ) − T (δT ∞ k ) is a positivedefinite 2-form with respect to δT * k .Consequently, the minimum of T is reached when δT k = δT ∞ k , and the system satisfies an incremental least action principle defined by the functional T .

C. An illustrative numerical example
We consider the following unidirectional and adimensional diffusion problem: with the initial and boundary conditions The solution of this problem can easily be obtained by separation of variables and is given by (the subscript "ana" meaning "analytic") with This problem has been solved by three different numerical methods.The first method is the finite-difference method using an implicit time scheme, and is denoted as "FD" below.
The second and third methods (IRM1 and IRM2) correspond to the incremental Ritz method (IRM) introduced at the end of Sec.III D and applied to the following functional: This functional is the analog of (54), with ∂ q = ∅, except that its minimum satisfies the energy balance rather than the entropy balance: It leads to a lighter implementation due to the absence of temperature at denominator, and is defined over the set of temperature increments δT (x) satisfying δT (0) = 0 and δT (1) = 0.The first implementation of the incremental Ritz method IRM1 has been carried out with a parabolic increment reminiscent of the initial condition: in which a 1 is a scalar parameter.The implementation of the second incremental Ritz method IRM2 is associated with the same kind of increment but enriched by its proper square; the temperature increment thus depends on two parameters b 1 and b 2 : For IRM1 and IRM2, the algorithm is started by initializing the heat flux q appearing in (63) by −∂ x T ini (x).For each time increment, the parameter a 1 or the couple (b 1 , b 2 ) are calculated by evaluating the functional (63) on the increments (66) and (67), and by solving to localize the minimum value of T (δT 1 ) and T (δT 2 ).At the end of each time increment, the heat flow is incremented by −∂ x δT • for the next evaluation of T (δT The resolution was made over 30 time increments δt = 0.01 starting from t = 0.00. To compare the performance of the different numerical methods, we counted the n degrees of freedom (DOFs) necessary to achieve the simulations since in practical situations, this parameter is known to have a great impact on the CPU time of a simulation.For the finite-difference method (FD), the number of DOFs is the number of nodes used to discretize the interval [0; 1] minus 2, because the temperature is given at the two extremal nodes x = 0 and x = 1.For the incremental Ritz method, the number of DOFs is 1 (IRM1) or 2 (IRM2).For each numerical method leading to a temperature field T (x), we introduce two error criteria measuring the discrepancy between T (x) and T ana (x) and between their derivatives with respect to x.The criteria compare T (x) and T ana (x) for five arbitrary time values (t 1 = 0.00, t 2 = 0.05, t 3 = 0.10, t 4 = 0.20, t 5 = 0.30): The two constants E 0 T and E 0 Q have been adjusted to have E T = E Q = 1 for n = 1.To compare the two implementations of the incremental Ritz method with the finite-difference method, we varied n and calculated the corresponding values of E T and E Q (see Fig. 2).On this example, a good performance of the two implementations of the incremental Ritz method is observed.Indeed, it is necessary to have approximately 6 and 10 times more DOFs for the FD method to have a comparable value of E T to the IRM method (about 10 and 50 times for E Q ).

V. CONCLUDING REMARKS
In this paper, we have shown that any thermodynamic potential satisfies a differential conservation law which, once generalized in the sense of extended irreversible thermodynamics and integrated over a finite volume, fulfills an incremental least action principle in time.In addition, a physical meaning of self-adjointness (a necessary and sufficient condition for the Lagrangian to exist) has been proposed as the symmetry of the Tisza and Onsager matrices (T and L).One of the implications of our work is that, using the incremental minimum principle presented in Sec.III D, it is possible to update the values of the coefficients in T and L when time is changed by an increment δt.Consequently, the minimum principle could be used in practice as an incremental Ritz method (IRM) to solve nonlinear problems for which T and L may depend on the extensive and intensive variables e and I as well as the generalized flows j and forces X .As shown with the linear problem of Sec.IV C, this method is efficient in considerably reducing the number of degrees of freedom associated with the numerical simulation.However, the formulation of a suitable increment δI to evaluate T [see Eqs. ( 29) and (36)] can be a challenging task.Indeed, δI has to satisfy boundary conditions, and can be much less intuitive to formulate in the presence of couplings between several different physical processes.
Most of the attempts to find a variational form of equations governing one or more dissipative processes classically superimpose an equilibrium term and a dissipative term.The equilibrium term is defined from a thermodynamic potential [27][28][29][30][31] and the dissipative term from integrals over time.The present approach is consistent with these approaches except that the formulation is incremental in time.This approach is closely related to the Lagrangian approach proposed by Moroz [29] in the framework of nonequilibrium chemical thermodynamics for which the Lagrangian is written as the sum of a thermodynamic potential and a dissipative Rayleigh function that is bilinear with respect to flows.Another very close link seems also to exist with the incremental variational principle proposed in [32].
An important application of the Lagrangian and Hamiltonian formalisms is the use of the Liouville theorem which states that the volume of the phase space associated with a set of massive points remains constant when time varies.For the particular case of a single oscillator of momentum p and position q, this volume is given by β = d p ∧ dq and is an absolute integral invariant (AII) for the generator ∂ t + q ∂ q + ṗ ∂ p .Such an AII has proven to be very useful in numerical simulations with the development of "symplectic numerical schemes."Indeed, these schemes are designed such that the invariant β is exactly preserved along numerical iterations, i.e., β ν+1 = β ν , with ν the index of time increment.Consequently, they exhibit stability properties that are significantly improved with respect to classical numerical schemes (such as implicit or explicit Euler, Runge-Kutta).Our approach proposed a new AII as a 5-form defined in Eq. ( 23).It will be of great interest to formulate numerical schemes for which is preserved in a numerical way, and to analyze their numerical stability.Work is currently in progress to achieve this goal.

L[X
The last equation being true ∀α, we conclude that Since the {e k } are independent 0-forms, we conclude that and finally, identifying terms in , obtain

APPENDIX B: MATHEMATICAL RESULTS
We recall here the notations used in this paper: ∧ stands for the exterior product, d for the exterior differential.If not otherwise stated, the Einstein summation convention on repeated indexes is adopted.We recall that for any smooth function ξ (a) of p variables {a k , k ∈ [[1; p]]}, we have and then d 2 ξ = 0 as soon as ξ satisfies the Maxwell relations on its partial derivatives.Furthermore, for any smooth function f of time t and space variables (denoted x, y, z throughout the paper), the following equality holds: The proof of that result is straightforward: We recall that is an absolute integral invariant for ζ if and only if ζ ⌟ and ζ ⌟ d vanish (see, e.g., [24]).We now calculate these two interior products.We have The only way to have ζ ⌟ = 0 is then to fulfill the condition T lk = T kl .For the second interior product, we first calculate the exterior differential of :