# Thermodynamically consistent Cahn-Hilliard-Navier-Stokes equations using the
metriplectic dynamics formalism

Azeddine Zaidni  Philip J Morrison  Saad Benjelloun

###### Abstract

Cahn-Hilliard-Navier-Stokes (CHNS) systems describes flows with two-phases,
e.g., a liquid with bubbles. Obtaining constitutive relations for general
dissipative processes for such a systems, which are thermodynamically
consistent, can be a challenge. We show how the metriplectic 4-bracket
formalism [30] achieves this in a straightforward, in fact algorithmic,
manner. First, from the noncanonical Hamiltonian formulation for the ideal
part of a CHNS system we obtain an appropriate Casimir to serve as the entropy
in the metriplectic formalism that describes the dissipation (e.g. viscosity,
heat conductivity and diffusion effects). General thermodynamics with the
thermodynamic conjugates of concentration and chemical potential are included.
Having expressions for the Hamiltonian (energy), entropy, and Poisson bracket,
we describe a procedure for obtaining a metriplectic 4-bracket that describes
thermodynamically consistent dissipative effects. The 4-bracket formalism
leads naturally to a general CHNS system that allows for anisotropic surface
energy effects. This general CHNS system reduces to cases in the literature,
to which we can compare.

###### keywords:

Two phase fluid flow , Cahn-Hilliard , thermodynamic consistency ,
metriplectic dynamics , Hamiltonian structure

###### PACS:

0000 , 1111

###### MSC:

0000 , 1111

††journal:

[inst1]organization=Mohammed VI Polytechnic University, College of Computing
lab,addressline=Lot 660, Hay Moulay Rachid , city=Ben Guerir, postcode=43150,
country=Morocco

[inst2]organization=Department of Physics and Institute for Fusion Studies,
University of Texas at Austin,city=Austin, postcode=78712, state=TX,
country=USA [inst3]organization=Makhbar Mathematical Sciences Research
Institute,city=Casablanca, country=Morocco

###### Contents

  1. 1 Introduction
  2. 2 Metriplectic framework and the Gibbs-Navier-Stokes system
    1. 2.1 Description of the Gibbs-Navier-Stokes system
    2. 2.2 Noncanonical Poisson bracket of the Gibbs-Euler system
    3. 2.3 Metriplectic 4-bracket for the Gibbs-Navier-Stokes system
      1. 2.3.1 General metriplectic 4-bracket dynamics
      2. 2.3.2 General Kulkarni-Nomizu construction
      3. 2.3.3 Metriplectic 4-bracket for the GNS system
    4. 2.4 GNS metriplectic 2-bracket and conventional fluxes and affinities
      1. 2.4.1 GNS metriplectic 2-bracket
      2. 2.4.2 Fluxes and affinities for the GNS system
  3. 3 The Cahn-Hilliard-Navier-Stokes system
    1. 3.1 Hamiltonian and Entropy functional forms
    2. 3.2 Noncanonical Poisson bracket of the Cahn-Hilliard-Euler system
    3. 3.3 Metriplectic 4-bracket for the Cahn-Hilliard-Navier-Stokes system
    4. 3.4 Metriplectic 2-bracket for the Cahn-Hilliard-Navier-Stokes system
  4. 4 Summary and Conclusions
  5. A Semidefinite curvature

##  1 Introduction

The well-known Navier-Stokes equations govern the motion of a single-phase
fluid. However, in the case of two-phase fluids, chemical reactions, changes
of phase, and migration between substances of phases become significant and
cannot be disregarded. J. W. Cahn and J. E. Hilliard were the first to
formulate the mathematical equations that describe phase separation in a such
a binary fluid [4]. Here we investigate generalizations that combine the Cahn-
Hilliard equation with equations that describe the dynamics of fluid flow,
referred to as Cahn-Hilliard-Navier-Stokes (CHNS) systems. CHNS systems aim to
describe the hydrodynamic properties of a mixture of two phases such as
bubbles in a liquid. To narrow down the already broad scope, we assume that
the two fluids share the same velocity field, yet we allow for both extended
thermodynamics and diffusive interfaces between the two phases.

A substantial hurdle in developing CHNS type systems, systems with a variety
of constitutive relations, is to ensure thermodynamic consistency, i.e.,
adherence to the first law of thermodynamics, which in this context is to
produce a set of dynamical equations that conserve energy, and the second law
which in this context means the dynamical production of entropy, ensuring the
relaxation (asymptotic stability) to thermodynamic equilibrium. Here we
propose an algorithm for constructing such systems, an algorithm that produces
a large set of CHNS systems.

The algorithm has four steps: i) Select a set of dynamical variables. For a
CHNS system these will be
$\psi\coloneqq\\{\mathbf{m}=\rho\mathbf{v},\rho,\tilde{c}=\rho c,\sigma=\rho
s\\}$, which are the momentum density, mass density, volume concentration of
one of the constituents, and entropy density, respectively. ii) The next step
is to select energy and entropy functionals, $H$ and $S$, dependent on the
dynamical variables. The choice of these functionals is based on the physics
of the phenomena one wishes to describe. iii) The third step of the algorithm
is to obtain the noncanonical Poisson bracket [see 25] of the ideal
(nondissipative) part of the theory that has the chosen entropy as a Casimir
invariant. Since the work of Morrison and Greene [28], Poisson brackets for a
great many systems, including fluid and magnetofluid systems, have been found
[e.g. 21, 36, 1, 7, 38]. Thus, this step may be immediate. Alternatively, it
may be achieved by a coordinate change from a known Hamiltonian theory in
order to align with the chosen entropy functional. In either case, we obtain
at this stage a noncanonical Hamiltonian system. iv) The final step is to
construct a metriplectic 4-bracket as described in Morrison and Updike [30].
Although there are standard metriplectic 4-bracket constructions, there is
freedom at this last step to describe a variety of types of dissipation.
However, a natural choice follows upon consideration of the form of an early
metriplectic bracket [23]. Given $H$, $S$, and the 4-bracket, the dynamical
system with thermodynamically consistent dissipation is produced.

We apply the algorithm to two cases. First, in § 2, we consider a system where
the fluid thermodynamics is extended by allowing the internal energy to depend
on a concentration variable, with the chemical potential being its
thermodynamic dual. Because Gibbs introduced the notion of chemical potential,
we refer to the Hamiltonian version of this fluid systems as the Gibbs-Euler
(GE) system and the dissipative version as the Gibbs-Navier-Stokes (GNS)
system. It is a thermodynamically consistent version of the compressible
Navier-Stokes equations with the inclusion of this concentration variable for
describing a second phase of the fluid. The GNS system generalizes the early
work of Eckart [8, 9] and the treatment in [16]; it allows for all possible
thermodynamics fluxes. Next, in § 3, a general form of CHNS system is
produced, a form that models surface tension effects and allows for diffuse
interfaces. Our work is motivated in large part by the substantial works of
Anderson et al. [2] and Guo and Lin [17], which we generalize by obtaining a
class of systems that includes theirs as special cases. There is a huge
literature on this topic and these papers contain many important references to
previous work. [Also, see 11, for a recent review.]

The GNS system of § 2 serves as a straightforward example of our algorithm. In
§§ 2.1 we describe the set of dynamical variables, properties of the system,
the energy and entropy functionals $H$ and $S$. This amounts to the first and
second steps of the algorithm. Then in §§ 2.2 the Hamiltonian formulation of
the dissipation free part of the system is presented. This is the third step
of the algorithm where the Poisson bracket is obtained, after a brief review
of the noncanonical Hamiltonian formalism. Given the early work of Morrison
and Greene [28] and the classification of extensions in [36], this step is
immediate. Based on the early and recent works of Morrison [23] and [30] the
fourth step of the algorithm is also immediate. In §§ 2.3 we first review the
metriplectic 4-bracket formalism and present the realization that applies for
the GNS system. Thus, the thermodynamically consistent GNS system is
determined. In §§ 2.4 we obtain the metriplectic 2-bracket equations of
motion, and the determined fluxes and affinities, making connection to
standard irreversible thermodynamics. Using the results of § 2, we proceed in
§ 3 to obtain the main result of the paper, our general CHNS system that can
describe diffuse interface effects. The first and second steps of our
algorithm are taken in §§ 3.1, while the third step, obtaining the correct
Poisson bracket, is undertaken in §§ 3.2. In order to complete this step, one
must find the Poisson bracket for which the entropy of the second step is a
Casimir invariant, which we find can be achieved by a simple coordinate
transformation. The fourth step of the algorithm is taken in §§ 3.3. Here a
choice of metriplectic 4-bracket gives a general class of thermodynamically
consistent CHNS systems, a class that contains previous results as special
cases. The formalism also shows how one can transform to a simple entropy
variable at the expense of a more complicated internal energy. As in § 2, in
§§ 3.4 we reduce to the metriplectic 2-bracket. Finally in § 4 we briefly
summarize and make a few comments about ongoing and future work.

##  2 Metriplectic framework and the Gibbs-Navier-Stokes system

In this section we describe general features of the metriplectic framework in
the context of the GNS system, a generalization of the Navier-Stokes equations
that includes the dual thermodynamical variables of concentration and chemical
potential.

###  2.1 Description of the Gibbs-Navier-Stokes system

The GNS for 2 phase flow proceeds on familiar ground [8, 9, 16]. It amounts to
the single phase thermodynamic Navier-Stokes system or as it is sometimes
called the Fourier Navier-Stokes system with the dispersed phase described by
the addition of a concentration variable, $c$, giving the set of dynamical
variables $\psi=\\{\mathbf{v},\rho,c,s\\}$. Here we review global aspects of
this known system, before showing how it emerges from the metriplectic
formalism.

We suppose the mixture of two phases are contained in a volume $\Omega$, and
we consider the following global quantities and their evolution:

$\displaystyle{M}$ $\displaystyle=\int_{\Omega}\rho\,,\ \ \dot{M}=0,$ (1)
$\displaystyle\mathbf{P}$ $\displaystyle=\int_{\Omega}\rho\,\mathbf{v}\,,\ \
\dot{\mathbf{P}}=-\int_{\partial\Omega}\bar{J}_{\mathbf{m}}\cdot\mathbf{n},$
(2) $\displaystyle H$
$\displaystyle=\\!\int_{\Omega}\frac{\rho}{2}\,|\mathbf{v}|^{2}+\rho\,u(\rho,s,c)\,,\
\ \dot{H}=\\!-\int_{\partial\Omega}{\mathbf{J}}_{e}\cdot\mathbf{n}\,,$ (3)
$\displaystyle{C}$ $\displaystyle=\\!\int_{\Omega}\rho\,c\,,\ \
\dot{C}=\\!-\int_{\partial\Omega}{\mathbf{J}}_{c}\cdot\mathbf{n},$ (4)
$\displaystyle{S}$ $\displaystyle=\int_{\Omega}\rho\,s\,,\ \
\dot{S}=-\int_{\partial\Omega}{\mathbf{J}}_{s}\cdot\mathbf{n}+\int_{\Omega}\dot{s}^{prod}\,.$
(5)

Here $\rho$ is the density of the mixture, $\mathbf{v}$ is the mass-averaged
velocity of the mixture, $s$ is the specific entropy, and the phase variable
$c$ is the specific concentration (dimensionless mass concentration) that
determines how much of the dispersed phase of the mixture is present at a
point $\mathbf{x}\in\Omega\subset\mathbb{R}^{3}$. The variable $\tilde{c}=\rho
c$ is the mass density of the dispersed phase. The functionals ${M}$, ${P}$,
$H$ and ${S}$ are the total mass, momentum, energy, and entropy of the
mixture, respectively, while ${C}$ is the total mass of one of the
constituents. For convenience we will omit the incremental volume element for
integrations over $\Omega$, i.e., $\int_{\Omega}=\int_{\Omega}d^{3}x$ and we
used an over dot to mean the total derivative $d/dt$. The local thermodynamics
of the mixture is described by $u(\rho,s,c)$, the internal energy per unit
mass. For convenience the gravitational force is not considered, although its
inclusion is straightforward.

Quantities in the time derivatives of the basic functionals are as follows:
$\mathbf{n}$ is the unit outward normal vector of the boundary
$\partial\Omega$, ${\mathbf{J}}_{c}$ is the phase field flux, which depends on
gradient of the chemical potential, ${\bar{J}}_{\mathbf{m}}$ is the stress
tensor – surface forces – due to pressure and viscosity, ${\mathbf{J}}_{e}$
the energy flux that contains the rate of work done by the surface forces
(external energy), the rate of heat transfer and the rate of diffusivity in
phase field (internal energy), ${\mathbf{J}}_{s}$ is the net entropy flux
through the boundary, and $\dot{s}^{prod}$ is the local rate of entropy
production. The second law of thermodynamics is expressed by the requirement
that $\dot{s}^{prod}$ is non-negative.

For the GNS system the fluxes are given by

$\displaystyle{\mathbf{J}}_{c}$ $\displaystyle=-\bar{D}\cdot\nabla\mu\,,$ (6)
$\displaystyle\bar{J}_{\mathbf{m}}$
$\displaystyle=p\,{\bar{I}}-\bar{\bar{\Lambda}}:\nabla\mathbf{v}\,,$ (7)
$\displaystyle\mathbf{J}_{e}$
$\displaystyle=-\mathbf{v}\cdot\bar{\bar{\Lambda}}:\nabla\mathbf{v}-\bar{\kappa}\cdot\nabla
T-\mu\bar{D}\cdot\nabla\mu\,,$ (8) $\displaystyle{\mathbf{J}}_{s}$
$\displaystyle=-\frac{\bar{\kappa}}{T}\cdot\nabla T\,,$ (9)
$\displaystyle\dot{s}^{prod}$
$\displaystyle=\frac{1}{T}\Big{[}\nabla\mathbf{v}:\bar{\bar{\Lambda}}:\nabla\mathbf{v}+\frac{1}{T}\,\nabla
T\cdot\bar{\kappa}\cdot\nabla T$ $\displaystyle\hskip
56.9055pt+\nabla\mu\cdot\bar{D}\cdot\nabla\mu\Big{]}\geq 0\,,$ (10)

where $p$ is the pressure, $T$ is the temperature, $\bar{I}$ is the unit
tensor, $\bar{\kappa}$ is the thermal conductivity tensor, $\bar{D}$ is the
diffusion tensor, which along with $\bar{\kappa}$ is assumed to be a symmetric
and positive definite 2-tensor, and $\mu$ is the chemical potential. We allow
the possibility that phenomenological quantities such as $\bar{\kappa}$ and
$\bar{D}$ can depend on the dynamical variables. Here, $\bar{\bar{\Lambda}}$
is the viscosity 4-tenor, the usual rank 4 isotropic Cartesian tensor given by

$\Lambda_{ijkl}=\eta\left(\delta_{il}\delta_{jk}+\delta_{jl}\delta_{ik}-\frac{2}{3}\delta_{ij}\delta_{kl}\right)+\zeta\,\delta_{ij}\delta_{kl}\,,$
(11)

with viscosity coefficients $\eta$ and $\zeta$ and $i,j,k$ and $l$ taking on
values 1,2,3. Note, we use boldface as in the fluxes
${\mathbf{J}}_{c},{\mathbf{J}}_{e}$, and ${\mathbf{J}}_{s}$ to denote vectors,
an over bar as in ${\bar{J}}_{\mathbf{m}}$ to denote rank-2 tensors, and a
double over bar as in $\bar{\bar{\Lambda}}$ to denote rank-4 tensors. A single
“$\,\cdot\,$" is used for neighboring contractions, e.g.,
$(\bar{D}\cdot\nabla\mu)_{i}=D_{ij}\partial_{j}\mu$ and we use the double dot
convention, e.g., for the stress tensor
$(\bar{\bar{\Lambda}}:\\!\nabla\mathbf{v})_{ij}=\Lambda_{ijkl}\partial_{k}v_{l}$,
where repeated indices are summed over.

The volume density variables are
$\psi=(\rho,\mathbf{m}\coloneqq\rho\mathbf{v},\sigma\coloneqq\rho
s,\tilde{c}\coloneqq\rho c)$, where $\mathbf{m}$ is the momentum density,
$\sigma$ is entropy per unit volume and $\tilde{c}$ is the concentration per
unit volume. The local energy per unit volume is given by

$e=\frac{\,|\mathbf{m}|^{2}}{2\rho}+\rho u(\rho,s,c)\,.$ (12)

From the specific internal energy, $u(\rho,s,c)$, we have the thermodynamic
relations

$\mathrm{d}u=T\mathrm{d}s+\frac{p}{\rho^{2}}\mathrm{d}\rho+\mu\mathrm{d}c\,,$
(13)

where

$T=\frac{\partial u}{\partial s}\,,\qquad p=\rho^{2}\frac{\partial
u}{\partial\rho}\,,\quad\quad\mu=\frac{\partial u}{\partial c}\,.$ (14)

Given the content of this section, we have established the first step of our
algorithm for the GNS system, the determination of the dynamical variables
$\psi=\\{\mathbf{m}=\rho\mathbf{v},\rho,\tilde{c}=\rho c,\sigma=\rho s\\}$ or
alternatively the set $(\mathbf{v},\rho,c,s)$, and the second step of our
algorithm by making the choices of Hamiltonian $H$ of (3) and entropy $S$ of
(5). In the next section, §§ (2.2), we proceed to the third step of the
algorithm by obtaining the Hamiltonian structure for this system. This system
without dissipation is the GE system.

###  2.2 Noncanonical Poisson bracket of the Gibbs-Euler system

Given that the mixture is assumed to be confined in the domain $\Omega$, the
Eulerian scalars (volume forms) $(\rho,\tilde{c},\sigma)$ are functions from
space-time $\Omega\mapsto\mathbb{R}\rightarrow\mathbb{R}$, while the vector
field $\mathbf{m}$ maps $\Omega\times\mathbb{R}\mapsto T\Omega$, where
$T\Omega$ stands for the tangent bundle of the manifold $\Omega$. We will
forgo formal geometric considerations and suppose our infinite-dimensional
phase space has coordinates $\psi=(\mathbf{m},\rho,\tilde{c},\sigma)$ and
observables are functionals that map $\psi\mapsto\mathbb{R}$ at each fixed
time. We will denote the space of such functionals by $\mathcal{B}$. Then a
Poisson bracket is a bilinear operator
$\mathcal{B}\times\mathcal{B}\mapsto\mathbb{R}$ that fulfills the Leibniz rule
and is a realization of a Lie algebra [see e.g. 33, chap. 14]. The Leibniz
rule follows from that for the variational or functional derivative of
$F\in\mathcal{B}$, defined by

$\delta F[\psi;\delta\psi)=\lim_{\epsilon\rightarrow
0}\frac{F(\psi+\epsilon\delta\psi)-F(\psi)}{\epsilon}=\left\langle\frac{\delta
F}{\delta\psi},\delta\psi\right\rangle\,,$

where $\langle\,,\,\rangle$ is a vector space pairing and $\delta
F/\delta\psi$ is the functional derivative. This expression can be viewed as
the directional derivative of a functional $F$ at $\psi$ in the direction
$\delta\psi$ [see, e.g., 25, for a formal review of these notions].

The appropriate Poisson bracket, defined on two functionals
$F,G\in\mathcal{B}$, for the GE system is the following:

$\displaystyle\\{F,G\\}$
$\displaystyle=-\int_{\Omega}\mathbf{m}\cdot\left[F_{\mathbf{m}}\cdot\nabla
G_{\mathbf{m}}-G_{\mathbf{m}}\cdot\nabla F_{\mathbf{m}}\right]$
$\displaystyle+\rho\left[F_{\mathbf{m}}\cdot\nabla
G_{\rho}-G_{\mathbf{m}}\cdot\nabla F_{\rho}\right]$
$\displaystyle+\sigma\left[F_{\mathbf{m}}\cdot\nabla
G_{\sigma}-G_{\mathbf{m}}\cdot\nabla F_{\sigma}\right]$
$\displaystyle+\tilde{c}\left[F_{\mathbf{m}}\cdot\nabla
G_{\tilde{c}}-G_{\mathbf{m}}\cdot\nabla F_{\tilde{c}}\right]\,,$ (15)

where we compactified our notation by defining
$F_{\mathbf{m}}\\!\\!\coloneqq\delta F/\delta\mathbf{m}$,
$F_{\rho}\\!\\!\coloneqq\delta F/\delta\rho$, etc., the functional derivatives
with respect to the various coordinates $\psi$. That this is the appropriate
Poisson bracket is immediate; it is the Lie-Poisson bracket originally given
by Morrison and Greene [28] with the addition of the last line of (15)
involving the concentration, another volume density variable $\tilde{c}$.
Adding such a dynamical variable is common place in the fluid modeling of
plasmas over the last decades and fits within the general theory for extension
given by Thiffeault and Morrison [36]. By construction we have a Poisson
bracket that is a bilinear, antisymmetric, and either by the extension theory
or a relatively easy direct calculation using the techniques of Morrison [21]
it can be shown to satisfy the Jacobi identity, i.e.,

$\\{\\{F,G\\},H\\}+\\{\\{H,F\\},G\\}+\\{\\{G,H\\},F\\}=0\,,$ (16)

for all $F,G,H\in\mathcal{B}$. The Leibniz property, which is required for the
Poisson bracket to generate a vector field, is built into the definition of
functional derivative.

Upon inserting any functional of $\psi$, say an observable $o$, into the
Poisson bracket its evolution is determined by

$\partial_{t}o=\\{o,{H}\\}\,,$ (17)

where the Hamiltonian functional is the total energy of the system, where we
rewrite (3) as follows:

${H}[\rho,\mathbf{m},\sigma,\tilde{c}]=\int_{\Omega}\\!e\,=\int_{\Omega}\frac{\,|\mathbf{m}|^{2}}{2\rho}+\rho
u\left(\rho,\frac{\sigma}{\rho},\frac{\tilde{c}}{\rho}\right).$ (18)

In (17) and henceforth we use the shorthand $\partial_{t}=\partial/\partial
t$. Using the following functional derivatives:

$\displaystyle H_{\rho}$
$\displaystyle=-{|\mathbf{v}|^{2}}/{2}+u+{p}/{\rho}-sT-c\mu\,,\quad
H_{\mathbf{m}}=\mathbf{v}\,,$ $\displaystyle H_{\sigma}$
$\displaystyle=T\,,\quad H_{\tilde{c}}=\mu\,,$ (19)

the bracket form of 17 gives the ideal two-phase flow system

$\displaystyle\partial_{t}\mathbf{v}$
$\displaystyle=\\{\mathbf{v},H\\}=-\mathbf{v}\cdot\nabla\mathbf{v}-\nabla
p/\rho\,,$ (20) $\displaystyle\partial_{t}\rho$
$\displaystyle=\\{\rho,H\\}=-\mathbf{v}\cdot\nabla\rho-\rho\,\nabla\cdot\mathbf{v}\,,$
(21) $\displaystyle\partial_{t}\tilde{c}$
$\displaystyle=\\{\tilde{c},H\\}=-\mathbf{v}\cdot\nabla\tilde{c}-\tilde{c}\,\nabla\cdot\mathbf{v}\,,$
(22) $\displaystyle\partial_{t}\sigma$
$\displaystyle=\\{\sigma,H\\}=-\mathbf{v}\cdot\nabla\sigma-\sigma\,\nabla\cdot\mathbf{v}\,.$
(23)

Here we have dropped surface terms arising from integration by parts and have
used
$\delta\rho(\mathbf{x})/\delta\rho(\mathbf{x}^{\prime})=\delta(\mathbf{x}-\mathbf{x}^{\prime})$.
Equations (20)–(22) can also be written easily using , e.g.,
${\mathrm{D}\rho}/{\mathrm{D}t}\coloneqq{\partial\rho}/{\partial
t}+\mathbf{v}\cdot\nabla\rho$. These equations comprise the GE system.

Casimir invariants are special functionals $\mathfrak{C}$ that satisfy

$\\{F,\mathfrak{C}\\}=0\qquad\forall F\in\mathcal{B}\,,$ (24)

and thus are constants of motion for any Hamiltonian. From (15) we obtain the
following equations that a Casimir functional $\mathfrak{C}$ must satisfy:

$\nabla\cdot(\rho\,\mathfrak{C}_{\mathbf{m}})=\nabla\cdot(\sigma\,\mathfrak{C}_{\mathbf{m}})=\nabla\cdot(\tilde{c}\,\mathfrak{C}_{\mathbf{m}})=0\,,$
(25)

and

$m_{j}\nabla\mathfrak{C}_{m_{j}}+\partial_{j}(\mathbf{m}\,\mathfrak{C}_{m_{j}})+\rho\,\nabla\mathfrak{C}_{\rho}+\sigma\,\nabla\mathfrak{C}_{\sigma}+\tilde{c}\,\nabla\mathfrak{C}_{\tilde{c}}=0\,,$
(26)

where we use the shorthand
$\delta\mathfrak{C}/\delta\mathbf{m}\coloneqq\mathfrak{C}_{\mathbf{m}}$,
$\delta\mathfrak{C}/\delta\rho\coloneqq\mathfrak{C}_{\rho}$, etc. and
summation of repeated indices is assumed. For the purpose at hand we assume
$\mathfrak{C}$ is independent of $\mathbf{m}$, yielding the single condition

$\rho\,\nabla\mathfrak{C}_{\rho}+\sigma\,\nabla\mathfrak{C}_{\sigma}+\tilde{c}\,\nabla\mathfrak{C}_{\tilde{c}}=0\,.$
(27)

Equation (27) is satisfied by

$\mathfrak{C}=\int_{\Omega}\mathcal{C}(\rho,\sigma,\tilde{c})$ (28)

for any $\mathcal{C}$ that is Euler homogeneous of degree one, i.e., satisfies

$\mathcal{C}(\lambda\rho,\lambda\sigma,\lambda\tilde{c})=\lambda\,\mathcal{C}(\rho,\sigma,\tilde{c})\,.$
(29)

The proof of this is straightforward.

To complete the third step of our algorithm, the entropy functional must be
chosen from the set of Casimir invariants. Writing the Euler homogeneous
integrand as

$\mathcal{C}(\rho,\sigma,\tilde{c})=\rho f(\sigma/\rho,\tilde{c}/\rho)$

it is clear that

$S=\int_{\Omega}\rho\,s=\int_{\Omega}\\!\sigma\,$ (30)

lies in our set of Casimirs. This quantity was first shown to be a Casimir for
the ideal fluid in [21] and used for the thermodynamically consistent Navier-
Stokes metriplectic system in [23]. We note in passing, for other theories
that might have a nontraditional dynamical equilibrium playing the role of
thermodynamic equilibrium, one may wish to choose another Casimir.

###  2.3 Metriplectic 4-bracket for the Gibbs-Navier-Stokes system

Now let us turn to our fourth and final step of the algorithm, construction of
the metriplectic 4-bracket. To this end we review the formalism of [30] in
general terms and then apply it to the GNS system, which is a generalization
of an example given in that work.

####  2.3.1 General metriplectic 4-bracket dynamics

The metriplectic 4-bracket theory was introduced by Morrison and Updike [30]
to describe the dissipative dynamic. Let us briefly recall the metriplectic
4-bracket description in infinite dimensions. In this description, we consider
the dynamics of classical field theories with multi-component fields

$\chi(z,t)=\left(\chi^{1}(z,t),\chi^{2}(z,t),\ldots,\chi^{M}(z,t)\right)$ (31)

defined on $z=(z^{1},z^{2},\dots,z^{N})$ for times $t\in\mathbb{R}$. Here we
use $z$ to be a label space coordinate with the volume element $d^{N}\\!z$,
but with the domain unspecified. In fluid mechanics this domain would be
$\Omega$, the 3-dimensional domain occupied by the fluid and recall we used
$\mathbf{x}$ for the coordinate of this point. In general we suppose that
$\chi^{1},\ldots,\chi^{M}$ are real-valued functions of $z$ and $t$. Given the
space of functionals of $\chi$, $\mathcal{B}$, we define 4-bracket as an
operator

$(\,\cdot\,,\,\cdot\,\,;\,\cdot\,,\,\cdot\,)\colon\mathcal{B}\times\mathcal{B}\times\mathcal{B}\times\mathcal{B}\rightarrow\mathcal{B}$
(32)

such that for any four functionals $F,K,G,N\in\mathcal{B}$ we have

$\displaystyle(F,G;K,N)$
$\displaystyle=\int\\!\\!d^{N}\\!z\\!\\!\int\\!\\!d^{N}\\!z^{\prime}\,\\!\\!\int\\!\\!d^{N}\\!z^{\prime\prime}\\!\\!\int\\!\\!d^{N}\\!z^{\prime\prime\prime}\,\hat{R}^{\alpha\beta\gamma\delta}$
$\displaystyle\times\frac{\delta F}{\delta\chi^{\alpha}(z)}\frac{\delta
G}{\delta\chi^{\beta}\left(z^{\prime}\right)}\frac{\delta
K}{\delta\chi^{\gamma}\left(z^{\prime\prime}\right)}\frac{\delta
N}{\delta\chi^{\delta}\left(z^{\prime\prime\prime}\right)}\,,$ (33)

where
$\hat{R}^{\alpha\beta\gamma\delta}(z,z^{\prime},z^{\prime\prime},z^{\prime\prime\prime})$
is a 4-tensor functional operator with coordinate form given by the following
integral kernel:

$\displaystyle\hat{R}^{\alpha\beta\gamma\delta}(z,z^{\prime},z^{\prime\prime},z^{\prime\prime\prime})[\chi]$
$\displaystyle=$ (34)
$\displaystyle\hat{R}(\mathbf{d}\chi^{\alpha}(z),\mathbf{d}\chi^{\beta}(z^{\prime}),\mathbf{d}\chi^{\gamma}(z^{\prime\prime}),\mathbf{d}\chi^{\delta}(z^{\prime\prime}))[\chi(z)]\,,$

where $\alpha,\beta,\gamma,\delta$ range over $1,2,\dots,M$. The 4-bracket is
assumed to satisfy the following proprieties:
(i) Linearity in all arguments, e.g, for all $\lambda\in\mathbb{R}$

$(F+\lambda H,K;G,N)=(F,K;G,N)+\lambda(H,K;G,N)$ (35)

(ii) The algebraic symmetries

$\displaystyle(F,K;G,N)$ $\displaystyle=-(K,F;G,N)$ (36)
$\displaystyle(F,K;G,N)$ $\displaystyle=-(F,K;N,G)$ (37)
$\displaystyle(F,K;G,N)$ $\displaystyle=(G,N;F,K)$ (38)

(iii) Derivation in all arguments, e.g.,

$(FH,K;G,N)=F(H,K;G,N)+(F,K;G,N)H\,.$ (39)

Here, as usual, $FH$ denotes point-wise multiplication. In addition, to ensure
entropy production we require

$\dot{S}=(S,H;S,H)\geq 0\,.$ (40)

Metriplectic 4-brackets that satisfy (35)–(40) are called minimal
metriplectic. In §§ 2.3.2 we will give a construction that ensures such
appropriate positive semidefiniteness.

The minimal metriplectic properties of metriplectic 4-brackets are reminiscent
of the algebraic properties possessed by a curvature tensor. In fact, every
Riemannian manifold naturally has a metriplectic 4-bracket, and $(S,H;S,H)$
provides a notion of sectional curvature [see 30].

From the metriplectic 4-bracket (33), the dissipative dynamics of an
observable $o$ is generated as follows:

$\displaystyle\partial_{t}o$
$\displaystyle=\left(o,H;S,H\right)=\int\\!\\!d^{N}\\!z\\!\\!\int\\!\\!d^{N}\\!z^{\prime}\,\\!\\!\int\\!\\!d^{N}\\!z^{\prime\prime}\\!\\!\int\\!\\!d^{N}\\!z^{\prime\prime\prime}\,\hat{R}^{\alpha\beta\gamma\delta}$
$\displaystyle\times\frac{\delta o}{\delta\chi^{\alpha}(z)}\frac{\delta
H}{\delta\chi^{\beta}(z^{\prime})}\frac{\delta
S}{\delta\chi^{\gamma}(z^{\prime\prime})}\frac{\delta
H}{\delta\chi^{\delta}(z^{\prime\prime\prime})}\,.$ (41)

If we choose $o$ to be the Hamiltonian $H$, then
$\dot{H}=\left(H,H;S,H\right)\equiv 0$ by the antisymmetry condition of (36).
If we choose $o$ to be the entropy $S$, then $\dot{S}=\left(S,H;S,H\right)\geq
0$ by (40).

The dissipative dynamics generated by 4-bracket on our set of field variables
$\chi$ is given by

$\displaystyle\partial_{t}{\chi}^{\alpha}(z)$
$\displaystyle=\left(\chi^{\alpha},H;S,H\right)$
$\displaystyle=\\!\\!\int\\!\\!d^{N}\\!z^{\prime\prime}\,G^{\alpha\beta}(z,z^{\prime\prime})\frac{\delta
S}{\delta\chi^{\beta}(z^{\prime\prime})}\,,$ (42)

where the $G$-metric is given as follows:

$\displaystyle G^{\alpha\gamma}(z,z^{\prime\prime})$
$\displaystyle\coloneqq\int\\!\\!d^{N}\\!z^{\prime}\,\\!\\!\int\\!\\!d^{N}\\!z^{\prime\prime}\,R^{\alpha\beta\gamma\delta}(z,z^{\prime},z^{\prime\prime},z^{\prime\prime\prime})$
$\displaystyle\hskip 28.45274pt\times\frac{\delta
H}{\delta\chi^{\beta}(z^{\prime})}\frac{\delta
H}{\delta\chi^{\delta}(z^{\prime\prime\prime})}\,.$ (43)

For the full metriplectic dynamics we would add the Poisson bracket
contribution to the above. Equation (42) is written so as to show that it
amounts to a gradient systems with the entropy $S$ as generator.

####  2.3.2 General Kulkarni-Nomizu construction

We can easily create specific metriplectic 4-brackets that have the minimal
metriplectic properties: the requisite symmetries and the positive
semidefiniteness $(S,H;S,H)$. We do this by using the Kulkarni-Nomizu (K-N)
product [20, 31]. See also Fiedler [13] for relevant theorems. Consistent with
the bracket formulation of (33), we deviate from the conventional K-N product
by working on the dual. Given two symmetric operator fields, say $\sum$ and
$M$, operating on the variational derivatives; we again use the subscript
notation when convenient,

$F_{\chi}\coloneqq\frac{\delta F}{\delta\chi}=\left(\frac{\delta
F}{\delta\chi^{1}},\frac{\delta F}{\delta\chi^{2}},\ldots,\frac{\delta
F}{\delta\chi^{M}}\right)\,,$

the K-N product is defined as follows:

$\displaystyle(\Sigma\wedge M)\left(dF,dK,dG,dN\right)$
$\displaystyle=\,\Sigma\left(dF,dG\right)M\left(dK,dN\right)$
$\displaystyle-\Sigma\left(dF,dN\right)M\left(dK,dG\right)$
$\displaystyle+M\left(dF,dG\right)\Sigma\left(dK,dN\right)$
$\displaystyle-M\left(dF,dN\right)\Sigma\left(dK,dG\right)\,.$ (44)

A finite-dimensional form of $\Sigma$ would be a symmetric contravariant
2-tensor, say $\gamma$, and this would give the term

$\gamma(df,dg)=\gamma^{ij}\frac{\partial f}{\partial z^{i}}\frac{\partial
g}{\partial z^{j}}\,.$ (45)

A conventional for K-N product would involve rank 2 covariant tensors. The
form of (45) suggests a general form in infinite dimensions would be

$\Sigma(dF,dG)=\\!\\!\int\\!\\!d^{N}\\!z\\!\\!\int\\!\\!d^{N}\\!z^{\prime}\,\Sigma^{\alpha\beta}(z,z^{\prime})\frac{\delta
F}{\delta\chi^{\alpha}(z)}\frac{\delta G}{\delta\chi^{\beta}(z^{\prime})}\,,$
(46)

where $\Sigma^{\alpha\beta}(z,z^{\prime})$ is symmetric in both $\alpha,\beta$
and $z,z^{\prime}$ and operates to the right on both functional derivatives.
For example,

$\Sigma^{\alpha\beta}(z,z^{\prime})=L_{ab}^{\alpha\beta}(z,z^{\prime})\mathcal{L}^{a}\mathcal{L}^{\prime
b}\,,$ (47)

where $L_{ab}^{\alpha\beta}$ is symmetric and $\mathcal{L}^{a}$ is a
differential operator. This implies, e.g.,

$\displaystyle\Sigma(dF,dG)$
$\displaystyle=\int\\!\\!d^{N}\\!z\\!\\!\int\\!\\!d^{N}\\!z^{\prime}\,L_{ab}^{\alpha\beta}(z,z^{\prime})$
$\displaystyle\hskip 28.45274pt\times\mathcal{L}^{a}\frac{\delta
F}{\delta\chi^{\alpha}(z)}\,\mathcal{L}^{\prime b}\frac{\delta
G}{\delta\chi^{\beta}(z^{\prime})}\,.$ (48)

With an expression for $M$ similar to (46), a term in the K-N decomposition
would have the following form:

$\displaystyle\int\\!\\!d^{N}\\!z\\!\\!\int\\!\\!d^{N}\\!z^{\prime}\,\\!\\!\int\\!\\!d^{N}\\!z^{\prime\prime}\\!\\!\int\\!\\!d^{N}\\!z^{\prime\prime\prime}\,\Sigma^{\alpha\beta}(z,z^{\prime})\,M^{\gamma\delta}(z^{\prime\prime},z^{\prime\prime\prime})$
$\displaystyle\qquad\qquad\times\ \frac{\delta
F}{\delta\chi^{\alpha}(z)}\frac{\delta
G}{\delta\chi^{\beta}(z^{\prime})}\,\frac{\delta
K}{\delta\chi^{\gamma}(z^{\prime\prime})}\frac{\delta
N}{\delta\chi^{\delta}(z^{\prime\prime\prime})}\,,$ (49)

which could be generalized further by adding filtering kernels.

It is easy to see that brackets constructed with this K-N product will have
all of the algebraic symmetries described in §§ 2.3.1. In addition, it is
shown in appendix A using the Cauchy-Schwarz inequality that positivity of
$(S,H;S,H)$ is satisfied, if both $\Sigma$ and $M$ are positive semidefinite.
Moreover, if one of $\Sigma$ or $M$ is positive definite, defining an inner
product, then the sectional curvature of (37) satisfies $(S,H;S,H)\geq 0$ with
equality if and only if $\delta S/\delta\chi\propto\delta H/\delta\chi$. Thus,
it is not difficult to build minimal metriplectic 4-brackets.

Alternative to (46) we can define $\Sigma(dF,dG)$ pointwise as

$\displaystyle\Sigma(dF,dG)(z)$
$\displaystyle\coloneqq\\!\int\\!\\!d^{N}\\!z^{\prime}\,\Sigma^{\alpha\beta}(z,z^{\prime})\frac{\delta
F}{\delta\chi^{\alpha}(z)}\frac{\delta G}{\delta\chi^{\beta}(z^{\prime})}$
$\displaystyle=A^{\alpha\beta}(z)\frac{\delta
F}{\delta\chi^{\alpha}(z)}\frac{\delta G}{\delta\chi^{\beta}(z)}\,,$ (50)

which could follow from (46) if we added an additional argument to $\Sigma$.
Then, with a corresponding form for $M$ the algebraic curvature symmetries
would be induced in the integrand. This is the case for our present purposes,
where we assume the specific K-N form given in [30], viz. where the 4-bracket
is given by

$(F,K;G,N)=\\!\\!\\!\int d^{N}\\!z\,W\,(\Sigma\wedge
M)\left(dF,dK,dG,dN\right)\,,$

where $W$ is an arbitrary weight, possibly depending on $\chi$ and $z$, that
multiplies $(\Sigma\wedge M)$ where all of the functional derivatives are
evaluated a the same point, $z$. (See (52) and (53) below.) In §§ 2.3.3 we
will see that this form of 4-bracket is sufficient for the CHNS system of
interest.

####  2.3.3 Metriplectic 4-bracket for the GNS system

Now suppose our multi-component field variable $\chi$ is that for the
multiphase fluid, i.e.,

$\psi(\mathbf{x},t)=(\mathbf{m}(\mathbf{x},t),\rho(\mathbf{x},t),\tilde{c}(\mathbf{x},t),{\sigma}(\mathbf{x},t))$
(51)

and consider a specific, but still quite general, form of the K-N
construction, one adaptable to the GNS type of system. For multi-component
fields $\psi$ of our fluid we could choose

$\displaystyle M(dF,dG)$
$\displaystyle=F_{\psi^{\gamma}}A^{\gamma\delta}G_{{\psi^{\delta}}}\,,$ (52)
$\displaystyle\Sigma(dF,dG)$ $\displaystyle=\nabla F_{\psi^{\alpha}}\cdot
B^{\alpha\beta}\cdot\nabla G_{{\psi^{\beta}}}\,,$ (53)

where repeated indices are to be summed, $A^{\gamma\delta}=A^{\delta\gamma}$,
and in coordinates $B_{ij}^{\alpha\beta}$ is symmetric in both $i,j=1,2,3$ and
$\alpha,\beta=1,\dots,6$, which are indices that range over the six fields of
$\psi$. Here, the nablas are contracted on $i$ and $j$. With the choices of
(52) and (53), the metriplectic 4-bracket is

$\displaystyle(F,K;G,N)$ $\displaystyle=\int_{\Omega}\nabla
F_{\psi^{\alpha}}\cdot B^{\alpha\beta}\cdot\nabla
G_{{\psi^{\beta}}}\,K_{\psi^{\gamma}}A^{\gamma\delta}N_{{\psi^{\delta}}}$
$\displaystyle-\nabla F_{\psi^{\alpha}}\cdot B^{\alpha\beta}\cdot\nabla
N_{{\psi^{\beta}}}\,K_{\psi^{\gamma}}A^{\gamma\delta}G_{{\psi^{\delta}}}$
$\displaystyle+\nabla K_{\psi^{\alpha}}\cdot B^{\alpha\beta}\cdot\nabla
N_{{\psi^{\beta}}}\,F_{\psi^{\gamma}}A^{\gamma\delta}G_{{\psi^{\delta}}}$
$\displaystyle-\nabla K_{\psi^{\alpha}}\cdot B^{\alpha\beta}\cdot\nabla
G_{{\psi^{\beta}}}\,F_{\psi^{\gamma}}A^{\gamma\delta}N_{{\psi^{\delta}}}\,.$
(54)

Entropy production is governed by

$\displaystyle(S,H;S,H)$ $\displaystyle=\int_{\Omega}\nabla
S_{\psi^{\alpha}}\cdot B^{\alpha\beta}\cdot\nabla
S_{{\psi^{\beta}}}\,H_{\psi^{\gamma}}A^{\gamma\delta}H_{{\psi^{\delta}}}$ (55)
$\displaystyle+\nabla H_{\psi^{\alpha}}\cdot B^{\alpha\beta}\cdot\nabla
H_{{\psi_{\beta}}}\,S_{\psi^{\gamma}}A^{\gamma\delta}S_{{\psi^{\delta}}}$
$\displaystyle-2\,\nabla H_{\psi^{\alpha}}\cdot B^{\alpha\beta}\cdot\nabla
S_{{\psi_{\beta}}}\,S_{\psi^{\gamma}}A^{\gamma\delta}H_{{\psi^{\delta}}}\,.$

From the general results of appendix A it follows that $(S,H;S,H)\geq 0$ for
this special case if $A^{\alpha\beta}$ and $B^{\alpha\beta}$ are positive
semidefinte.

Observe that (53) could be replaced by the more general expression

$\Sigma(dF,dG)=\mathcal{L}F_{\psi^{\alpha}}\cdot
B^{\alpha\beta}\cdot\mathcal{L}G_{{\psi^{\beta}}}\,,$ (56)

where $\mathcal{L}$ is contained within a general class of pseudodifferential
operators. Later we will see an example of this.

Now consider an even more restrictive K-N product, a special case of (73) with
what appears to be the simplest K-N options. As discussed earlier, we do not
expect our 4-bracket to depend on functional derivatives with respect to
$\rho$, which could produce density diffusion. Thus, for $M$ we take

$M(dF,dG)=F_{{\sigma}}G_{{\sigma}}.$ (57)

The placement of the $\nabla$ in (73) leads to a diffusive type of relaxation,
so this is natural, and the simplest case would be to select $\Sigma$ with no
cross terms, i.e.,

$\displaystyle{\Sigma}(dF,dG)$ $\displaystyle=\nabla
F_{{\mathbf{m}}}:\bar{\bar{\Lambda}}_{1}:\nabla G_{{\mathbf{m}}}+\nabla
F_{{\sigma}}\cdot\bar{\Lambda}_{2}\cdot\nabla G_{{\sigma}}$
$\displaystyle\hskip
56.9055pt+\nabla{F}_{\tilde{c}}\cdot\bar{\Lambda}_{3}\cdot\nabla{G}_{\tilde{c}}\,,$
(58)

where the 4-tensor $\bar{\bar{\Lambda}}_{1}$ and the symmetric 2-tensors
$\bar{\Lambda}_{2}$ and $\bar{\Lambda}_{3}$ are to be determined. We make the
following choices

$\bar{\bar{\Lambda}}_{1}=\frac{\bar{\bar{\Lambda}}}{T},\qquad\bar{\Lambda}_{2}=\frac{\bar{\kappa}}{T^{2}},\qquad\bar{\Lambda}_{3}=\frac{\bar{D}}{T}\,,$
(59)

where $\bar{\bar{\Lambda}}$ is the isotropic Cartesian 4-tensor given by (11),
$\bar{\Lambda}_{2,3}$ are symmetric positive definite 2-tensors defined by the
previously introduced $\bar{\kappa}$ and $\bar{D}$. We take the weight $W$ to
be the Lagrange multiplier defined in §§ 2.4.1 i.e. $W=1$. Then, the 4-bracket
reads

$\displaystyle(F,K;G,N)$ $\displaystyle=$ (60)
$\displaystyle\int_{\Omega}\\!\frac{1}{T}\Big{[}\left[K_{\sigma}\nabla
F_{\mathbf{m}}-F_{\sigma}\nabla
K_{\mathbf{m}}\right]:\bar{\bar{\Lambda}}:\left[N_{\sigma}\nabla
G_{\mathbf{m}}-G_{\sigma}\nabla N_{\mathbf{m}}\right]$ $\displaystyle+\
\frac{1}{T}\,\big{[}K_{\sigma}\nabla F_{{\sigma}}-F_{\sigma}\nabla
K_{{\sigma}}\big{]}\cdot\bar{\kappa}\cdot\big{[}N_{\sigma}\nabla
G_{{\sigma}}-G_{\sigma}\nabla N_{{\sigma}}\big{]}$
$\displaystyle+\big{[}K_{\sigma}\nabla F_{\tilde{c}}-F_{\sigma}\nabla
K_{\tilde{c}}\big{]}\cdot\bar{D}\cdot\big{[}N_{\sigma}\nabla
G_{\tilde{c}}-G_{\sigma}\nabla N_{\tilde{c}}\big{]}\Big{]}\,.$

Upon insertion of $H$ as given by (18) and $S$ from the set of Casimirs of §§
2.2 to be as in (30), the dynamics is given by

$\partial_{t}\psi^{\alpha}=\\{\psi^{\alpha},H\\}+(\psi^{\alpha},H;S,H)\,.$
(61)

Using ${H}_{\mathbf{m}}=\mathbf{v}$, ${H}_{\sigma}=T$, and ${S}_{\sigma}=1$,
the following GNS system is produced:

$\displaystyle\partial_{t}\mathbf{v}$
$\displaystyle=\\{\mathbf{v},H\\}+(\mathbf{v},H;S,H)$
$\displaystyle=-\mathbf{v}\cdot\nabla\mathbf{v}-\nabla
p/\rho+\frac{1}{\rho}\nabla\cdot(\bar{\bar{\Lambda}}:\nabla\mathbf{v})\,,$
(62) $\displaystyle\partial_{t}\rho$ $\displaystyle=\\{\rho,H\\}+(\rho,H;S,H)$
$\displaystyle=-\mathbf{v}\cdot\nabla\rho-\rho\,\nabla\cdot\mathbf{v}\,,$ (63)
$\displaystyle\partial_{t}\tilde{c}$
$\displaystyle=\\{\tilde{c},H\\}+(\tilde{c},H;S,H)$
$\displaystyle=-\mathbf{v}\cdot\nabla\tilde{c}-\tilde{c}\,\nabla\cdot\mathbf{v}+\nabla\cdot(\bar{D}\cdot\nabla\mu)\,,$
(64) $\displaystyle\partial_{t}\sigma$
$\displaystyle=\\{\sigma,H\\}+(\sigma,H;S,H)$
$\displaystyle=-\mathbf{v}\cdot\nabla\sigma-\sigma\,\nabla\cdot\mathbf{v}$
$\displaystyle\qquad+\nabla\cdot\left(\frac{\bar{\kappa}}{T}\cdot\nabla
T\right)+\frac{1}{T^{2}}\nabla T\cdot\bar{\kappa}\cdot\nabla T$
$\displaystyle\qquad+\frac{1}{T}\nabla\mathbf{v}:\bar{\bar{\Lambda}}:\nabla\mathbf{v}+\frac{1}{T}\nabla\mu\cdot\bar{D}\cdot\nabla\mu\,.$
(65)

By construction we automatically have energy conservation, i.e., for (18)
$\dot{H}=0$, and entropy production

$\displaystyle\dot{{S}}$ $\displaystyle=({S},{H};{S},{H})$
$\displaystyle=\int_{\Omega}\frac{1}{T}\bigg{[}\nabla\mathbf{v}:\bar{\bar{\Lambda}}:\nabla\mathbf{v}+\frac{1}{T}\nabla
T\cdot\bar{\kappa}\cdot\nabla T$ $\displaystyle\hskip
88.2037pt+\nabla\mu\cdot\bar{D}\cdot\nabla\mu\bigg{]}\geq 0\,.$ (66)

###  2.4 GNS metriplectic 2-bracket and conventional fluxes and affinities

For completeness we demonstrate two things in this subsection: how the
metriplectic 4-bracket formalism relates to the original binary metriplectic
formalism given in Morrison [22, 23, 24] and how it relates to conventional
nonequilibrium thermodynamics, making the connection between the 4-bracket K-N
construction and the phenomenology of thermodynamics fluxes and affinities
(sometimes called thermodynamic forces).

####  2.4.1 GNS metriplectic 2-bracket

As noted above we are concerned with the metriplectic dynamics introduced in
Morrison [22, 23, 24] [see also 26, 6], but we mention that other binary
brackets for describing dissipation were presented over the years [e.g. 19,
18, 29, 14, 15, 3, 10]. In addition we mention a recent alternative approach
to multiphase fluids, one based on constrained variational principles, that is
given in Eldred and Gay-Balmaz [12]. We refer the reader to Morrison and
Updike [30] for comparisons with other formulations and how they emerge from
the metriplectic 4-bracket.

Metriplectic dynamics was introduced as a means of building thermodynamically
consistent theories in terms of a binary bracket, which we now call the
metriplectic 2-bracket. The theory applies to a wide class of dynamical
systems, including both ordinary and partial differential equations. Evolution
of an observable $o$ using the metriplectic 2-bracket has the following form:

$\displaystyle\partial_{t}o=\\{o,\mathcal{F}\\}-(o,\mathcal{F})_{H}\,,$ (67)

where as before $\\{\,.\,\\}$ is the noncanonical Poisson bracket that
generates the ideal part of the dynamics, while now $(F,G)_{H}$, the
metriplectic 2-bracket, generates the dissipative part. The functional
$\mathcal{F}$ represents the global Helmholtz free energy of the system, and
is given by:

$\displaystyle\mathcal{F}={H}-\mathcal{T}{S}\,,$ (68)

where again ${H}$ is the Hamiltonian and ${S}$ the entropy selected from the
set of Casimirs of the noncanonical Poisson bracket (ensuring
$\\{F,{S}\\}=0\,$ for any functional $F$), and $\mathcal{T}$ is a uniform
nonnegative constant (a global temperature). The metriplectic 2-bracket
$(\,,\,)$ is assumed to be bilinear, symmetric, and satisfies

$\displaystyle(F,{H})_{H}\equiv 0\quad\text{ for any functional }F.$ (69)

Thus, metriplectic systems are thermodynamically consistent:

First law

(energy conservation):

$\displaystyle\dot{H}$
$\displaystyle=\\{{H},\mathcal{F}\\}-({H},\mathcal{F})_{H}$
$\displaystyle=\\{{H},{H}\\}+\mathcal{T}({H},{S})_{H}=0\,;$ (70)

Second law

(entropy production):

$\displaystyle\dot{S}$ $\displaystyle=\\{S,\mathcal{F}\\}-(S,\mathcal{F})_{H}$
$\displaystyle=-({S},{H})_{H}+\mathcal{T}({S},{S})_{H}=\mathcal{T}({S},{S})_{H}\geq
0\,,$ (71)

which follows because $\\{S,\mathcal{F}\\}\equiv 0$ and $({S},{H})_{H}\equiv
0$. As shown in Morrison and Updike [30] the metriplectic 2-bracket emerges
from the 4-bracket as follows:

$(F,{G})_{{H}}=(F,{H};{G},{H})\,,$ (72)

where for convenience here and henceforth we set $\mathcal{T}=1$. Because of
the minimal metriplectic properties of the 4-bracket, we are assured to have
the thermodynamic consistency of (70) and (71).

The 2-bracket that emerges from the general 4-bracket of (54) is the
following:

$\displaystyle(F,G)_{H}$ $\displaystyle=(F,H;G,H)$
$\displaystyle=\int_{\Omega}\nabla F_{\psi^{\alpha}}\cdot
B^{\alpha\beta}\cdot\nabla
G_{{\psi^{\beta}}}\,H_{\psi^{\gamma}}A^{\gamma\delta}H_{{\psi^{\delta}}}$
$\displaystyle-\nabla F_{\psi^{\alpha}}\cdot B^{\alpha\beta}\cdot\nabla
H_{{\psi^{\beta}}}\,H_{\psi^{\gamma}}A^{\gamma\delta}G_{{\psi^{\delta}}}$
$\displaystyle+\nabla H_{\psi^{\alpha}}\cdot B^{\alpha\beta}\cdot\nabla
H_{{\psi_{\beta}}}\,F_{\psi^{\gamma}}A^{\gamma\delta}G_{{\psi^{\delta}}}$
$\displaystyle-\nabla H_{\psi^{\alpha}}\cdot B^{\alpha\beta}\cdot\nabla
G_{{\psi_{\beta}}}\,F_{\psi^{\gamma}}A^{\gamma\delta}H_{{\psi^{\delta}}}\,,$
(73)

which in light of the K-N product satisfies $(F,H)_{H}=0$ for all $F$. This
will be true for any choice of the Hamiltonian $H$. Indeed, a special case of
this was used in (60) to obtain the thermodynamically consistent set of
equations (62)-(64). Recall, for this case $M$ and $\Sigma$ were chosen as in
(57) and (58), special cases of (52) and (53).

Another 2-bracket can be obtained from (54) by making a convenient choice of
variables; viz., instead of the variables $\psi$ of (134) in (73) we choose
the following:

$\xi(\mathbf{x},t)\coloneqq(\mathbf{m}(\mathbf{x},t),\rho(\mathbf{x},t),\tilde{c}(\mathbf{x},t),e(\mathbf{x},t))\,,$
(74)

where the total energy density is used instead of the entropy density as one
of our dynamical variables. That this is possible is well known in
thermodynamics because the entropy must be a monotonic increasing function of
the internal energy, which allows via the inverse function theorem
transformation between the extensive energy or extensive entropy
representations [5]. We will denote these density variables in order by
$\xi_{\alpha}$, $\alpha=1,\dots,6$. With this choice the Hamiltonian (18) is
given by

$H=\int_{\Omega}e=\int_{\Omega}\xi_{6}\,,$ (75)

$\delta H/\delta\xi_{\alpha}=\delta_{\alpha 6}$, the Kronecker delta, and
$\nabla\delta H/\delta\xi_{\alpha}\equiv 0$. Thus, (73) reduces to

$(F,G)_{H}=(F,H;G,H)=\int_{\Omega}\nabla F_{\xi_{\alpha}}\cdot
L_{\alpha\beta}\cdot\nabla G_{{\xi_{\beta}}}\,,$ (76)

where without loss of generality we set $A^{66}=1$ and
$B^{\alpha\beta}=L_{\alpha\beta}$ for this case. In §§ 2.4.2 we will follow
Coquinot and Morrison [6] and show how the bracket of (76) fits into the
framework of conventional nonequilibrium thermodynamics as, e.g., described in
de Groot and Mazur [16]. In this way we will physically identify the meaning
of $L_{\alpha\beta}$.

####  2.4.2 Fluxes and affinities for the GNS system

A fundamental equations of nonequilibrium thermodynamics is the general
thermodynamic identity

$\text{d}\sigma=X^{\alpha}\text{d}\xi_{\alpha},$ (77)

which relates $\sigma$, the entropy density, to the $\xi_{\alpha}$ densities
associated with conserved extensive properties and to
$X^{\alpha}\coloneqq\partial\sigma/\partial\xi_{\alpha}$, quantities called
affinities (or thermodynamic forces). All the densities are characterized by
the following conservation equations:

$\partial_{t}\xi_{\alpha}+\nabla\cdot\mathcal{J}_{\alpha}=0\,,$ (78)

where $\mathcal{J}_{\alpha}$ is at present an unknown flux associated with the
density $\xi_{\alpha}$. Then, the evolution of the entropy is given by

$\partial_{t}\sigma+\nabla\cdot\left(X^{\alpha}\mathcal{J}_{\alpha}\right)=\mathcal{J}_{\alpha}\cdot\nabla
X^{\alpha}\,.$ (79)

The righthand side of (79) is the dissipative term, which is the sum of the
fluxes $\mathcal{J}_{\alpha}$ contracted with $\nabla X^{\alpha}$. The linear
assumption of nonequilibrium processes amounts to relating fluxes and
affinities according to

$\mathcal{J}_{\alpha}=L_{\alpha\beta}X^{\beta}\,.$ (80)

If we identify the $L_{\alpha\beta}$ of (80) with that of (76), we see how
metriplectic brackets are related to the flux-affinity relations. Onsager
symmetry, assumed to arise from microscopic reversibility, amounts to the
symmetry $L_{\alpha\beta}=L_{\beta\alpha}$ and the semi-definiteness property
assures the second law, i.e., entropy growth. It remains to identify how the
fluxes of §§ 2.1 enter the picture.

To further identify the meaning of $L_{\alpha\beta}$ we revisit the
thermodynamics of (13), in light of our choice of the variables $\xi$ of
(134). Thus, we rewrite the thermodynamic relation 13 upon changing variables,

$T\text{d}\sigma=\text{d}e-\mathbf{v}\cdot\text{d}\mathbf{m}-\mu\text{d}\tilde{c}-g\text{d}\rho\,,$
(81)

where $e$ is the energy density of (12) and $g$ is a modified specific Gibbs
free energy, viz.

$g:=u-Ts+{p}/{\rho}-\mu c-{|\mathbf{v}|^{2}}/{2}\,.$ (82)

We have assumed in (1) that there is no flux associated with $\rho$, i.e., in
CHNS chemical reactions and/or particle creation and annihilation are ignored.
Thus, the phase space for the thermodynamics is smaller than that for the
Hamiltonian dynamics, because the variable $\rho$ as seen e.g. in (63) has no
dissipative terms. This leads us to focus on the thermodynamic variables
$(e,\mathbf{m},\tilde{c})$ and (81) reduces to

$T\mathrm{d}\sigma=\mathrm{d}e-\mathbf{v}\cdot\mathrm{d}\mathbf{m}-\mu\mathrm{d}\tilde{c}\,.$
(83)

Comparison of (77) and (83) suggests we require the affinities associated with
$\mathbf{m}$, $e$, and $\tilde{c}$. The conventional choices for these
affinities are $\nabla(1/T)$, $\nabla(-\mathbf{v}/T)$, and $\nabla(-\mu/T)$,
respectively [16]. However, examination of (10) or (66) suggests using instead
$\nabla T$, $\nabla\mathbf{v}$, and $\nabla\mu$, as was done in Coquinot and
Morrison [6].

The relationship between the flux-affinity relations in terms of these two
choices of bases are given by the following:

$\displaystyle\bar{J}_{\mathbf{m}}$
$\displaystyle=L_{\mathbf{m}e}\cdot\nabla\left({\frac{1}{T}}\right)+L_{\mathbf{m}\mathbf{m}}:\nabla\left(\frac{-\mathbf{v}}{T}\right)+L_{\mathbf{m}\tilde{c}}\cdot\nabla\left({\frac{-\mu}{T}}\right)$
$\displaystyle\hskip 14.22636pt=-\bar{\bar{\Lambda}}:\nabla\mathbf{v}\,,$ (84)
$\displaystyle{\mathbf{J}}_{e}$
$\displaystyle=L_{ee}\cdot\nabla\left(\frac{1}{T}\right)+L_{e\mathbf{m}}:\nabla\left(\frac{-\mathbf{v}}{T}\right)+L_{e\tilde{c}}\cdot\nabla\left(\frac{-\mu}{T}\right)$
$\displaystyle\hskip
14.22636pt=-\mathbf{v}\cdot(\bar{\bar{\Lambda}}:\nabla\mathbf{v})-\bar{\kappa}\cdot\nabla
T-\mu\bar{D}\cdot\nabla\mu\,,$ (85) $\displaystyle{\mathbf{J}}_{\tilde{c}}$
$\displaystyle=L_{\tilde{c}{e}}\cdot\nabla\left(\frac{1}{T}\right)+L_{\tilde{c}\mathbf{m}}:\nabla\left(\frac{-\mathbf{v}}{T}\right)+L_{\tilde{c}\tilde{c}}\cdot\nabla\left(\frac{-\mu}{T}\right)$
$\displaystyle\hskip 14.22636pt=-\bar{D}\cdot\nabla\mu\,.$ (86)

Recall $\bar{J}_{\mathbf{m}}$ is a 2-tensor, thus
$L_{\mathbf{m}e}=L_{e\mathbf{m}}$ is a 3-tensor, $L_{\mathbf{m}\mathbf{m}}$ is
a 4-tensor, and $L_{\mathbf{m}\tilde{c}}=L_{\tilde{c}\mathbf{m}}$ is a
3-tensor. Since $\mathbf{J}_{e}$ and $\mathbf{J}_{\tilde{c}}$ are vectors,
$L_{ee},L_{{e}\tilde{c}}=L_{\tilde{c}e}$, and $L_{\tilde{c}\tilde{c}}$ are
2-tensors. From (84), (85), and (86), we identify the components of
$L_{\alpha\beta}$ as follows:

$\displaystyle L_{\mathbf{m}e}=T\bar{\bar{\Lambda}}\cdot\mathbf{v}\,,\ \ \
L_{\mathbf{m}\mathbf{m}}=T\bar{\bar{\Lambda}}\,,\ \ \
L_{\mathbf{m}\tilde{c}}=0\,,\ \ \ L_{{e}\tilde{c}}=T\mu\bar{D}^{\prime}\,,$
$\displaystyle
L_{ee}=T^{2}\bar{\kappa}+T\mathbf{v}\cdot\bar{\bar{\Lambda}}\cdot\mathbf{v}+T\mu^{2}\bar{D}\,,\quad
L_{\tilde{c}\tilde{c}}=T\bar{D}\,.$ (87)

The metriplectic 2-bracket in terms of the $\xi$-variables is given by

$\displaystyle(F,G)_{H}$ $\displaystyle=\int_{\Omega}T\,\Big{[}\nabla
F_{e}\cdot(T\bar{\kappa}+\mathbf{v}\cdot\bar{\bar{\Lambda}}\cdot\mathbf{v}+\mu^{2}\bar{D})\cdot\nabla
G_{e}$ $\displaystyle+\nabla F_{\mathbf{m}}:\bar{\bar{\Lambda}}:\nabla
G_{\mathbf{m}}+\nabla F_{\tilde{c}}\cdot\bar{D}\cdot\nabla G_{\tilde{c}}$
$\displaystyle+\nabla F_{e}\cdot(\bar{\bar{\Lambda}}\cdot\mathbf{v}):\nabla
G_{\mathbf{m}}+\nabla G_{e}\cdot(\bar{\bar{\Lambda}}\cdot\mathbf{v}):\nabla
F_{\mathbf{m}}$ $\displaystyle+\mu\,\nabla F_{e}\cdot\bar{D}\cdot\nabla
G_{\tilde{c}}+\mu\,\nabla G_{e}\cdot\bar{D}\cdot\nabla
F_{\tilde{c}}\Big{]}\,.$ (88)

Upon writing

$S=\int_{\Omega}\sigma(\rho,e,\tilde{c},\mathbf{m})$ (89)

and using standard thermodynamic manipulations we obtain

$\displaystyle S_{e}={1}/{T}\,,\ \ S_{\mathbf{m}}=-{\mathbf{v}}/{T}\,,\ \
\mathrm{and}\ \ {S}_{\tilde{c}}=-{\mu}/{T}\,.$ (90)

Inserting these into $(o,S)_{H}$ using the 2-bracket of 88 yields the
dissipative terms of (62), (64), and with the manipulations of transforming
from $e$ to $\sigma$, those of the entropy equation (65). By direct
calculation, as well as by construction, we obtain $(H,S)_{H}=0$ and
$\dot{S}=({S},{S})_{H}=(S,H;S,H)\geq 0$ which reproduces (66).

To close the circle we transform the bracket of 88 in terms of the variables
$(e,\mathbf{m},\tilde{c})$ to one in terms of $(\sigma,\mathbf{m},\tilde{c})$
via the following chain rule formulas:

$\displaystyle G_{e}$ $\displaystyle\rightarrow G_{\sigma}/T,\qquad
G_{\mathbf{m}}\rightarrow G_{\mathbf{m}}-\mathbf{v}\,,G_{\sigma}/{T}\,,$
$\displaystyle\quad G_{\tilde{c}}\rightarrow
G_{\tilde{c}}-{\mu}\,G_{\sigma}/T\,.$ (91)

This calculation gives precisely the bracket of (60).

##  3 The Cahn-Hilliard-Navier-Stokes system

Now we apply our algorithm to obtain the Cahn-Hilliard-Navier-Stokes system
(CHNS) which allows for diffuse-interfaces. We follow the steps in order, just
as in § 2. However, here we have the additional step of aligning the desired
entropy functional with the Poisson bracket, so that it is indeed a Casimir
invariant.

###  3.1 Hamiltonian and Entropy functional forms

The phenomenon of material transport along an interface is known as the
Marangoni effect. The presence of a surface tension gradient naturally induces
the migration of particles, moving from regions of low tension to those of
high tension. This gradient can be triggered by a concentration gradient (or
also a temperature gradient). In two-phase theory the interface between phases
is regarded as being diffuse. According to the work of Taylor and Cahn [35],
one can model the diffuse interface by a single order parameter, say $\phi$
and with a free energy functional,

$\mathfrak{F}=\int_{\Omega}\frac{\epsilon}{2}\,\Gamma^{2}(\nabla\phi)+\frac{1}{\epsilon}V(\phi)\,,$
(92)

with $\Gamma$ being a homogeneous function of degree one, further details on
this will be provided later. Here $V$ can be any non-negative function that
equals zero at $\phi=\pm 1$ and $\epsilon$ is a small parameter that goes to
zero in the sharp-interface limit. We choose the order parameter $\phi$ to be
the concentration.

In the isotropic surface energy case Guo and Lin [17] develop a phase-field
model for two-phase flow, which is thermodynamically consistent. The modeling
based on a non-classical choice of energy and entropy, given respectively by

$\displaystyle H_{GL}$
$\displaystyle=\int_{\Omega}\frac{\rho}{2}\,|\mathbf{v}|^{2}+\rho
u(\rho,s,c)+\frac{\rho}{2}\lambda_{u}|\nabla c|^{2}\,,$ (93) $\displaystyle
S_{GL}$ $\displaystyle=\int_{\Omega}\rho s+\frac{\rho}{2}\lambda_{s}|\nabla
c|^{2}\,,$ (94)

where $u$ and $s$ stand for the classical specific internal energy and
entropy, respectively, while the coefficients $\lambda_{s}$ and $\lambda_{u}$
are constant parameters.

Alternatively, Anderson et al. [2] propose a model of phase-field of
solidification with convection, the model permits the interface to have an
anisotropic surface energy. The choice of energy and entropy are given by

$\displaystyle H_{AMW}$
$\displaystyle=\int_{\Omega}\frac{\rho}{2}\,|\mathbf{v}|^{2}+\rho
u(\rho,s,c)+\frac{\epsilon_{E}^{2}}{2}\Gamma^{2}(\nabla c)\,,$ (95)
$\displaystyle S_{AMW}$ $\displaystyle=\int_{\Omega}\rho
s-\frac{\epsilon_{S}^{2}}{2}\Gamma^{2}(\nabla c)\,,$ (96)

where the coefficients $\epsilon_{S}$ and $\epsilon_{E}$ are assumed to be
constant and $\Gamma$ is a homogeneous function of degree one that takes a
vector to a real number.

In this section, we explore a choice of energy and entropy functionals, from
which the previously mentioned choices are special cases, and we consider the
associated free energy functional, viz.

$\displaystyle H^{a}$
$\displaystyle=\int_{\Omega}\frac{\rho}{2}|\mathbf{v}|^{2}+\rho
u+\frac{\rho^{a}}{2}\lambda_{u}\Gamma^{2}(\nabla
c)\eqqcolon\\!\int_{\Omega}e_{\mathrm{Total}}^{a}\,,$ (97)
$\displaystyle{S^{a}}$ $\displaystyle=\int_{\Omega}\rho
s+\frac{\rho^{a}}{2}\lambda_{s}\Gamma^{2}(\nabla
c)\eqqcolon\\!\int_{\Omega}\sigma_{\mathrm{Total}}^{a}\,,$ (98)
$\displaystyle\mathfrak{F}^{a}\ $ $\displaystyle=\int_{\Omega}\rho
f+\frac{\rho^{a}}{2}\lambda_{f}(T)\Gamma^{2}(\nabla c)\,,$ (99)

where $u$, $s$ and $f$ stand for the classical specific internal energy,
entropy, and free energy, respectively, the coefficients $\lambda_{s}$ and
$\lambda_{u}$ are constant parameters, and $\lambda_{f}(T)$ is a parameter
depending on the temperature that will lead to anisotropic surface energy
effects. He have defined the total densities $e^{a}_{\mathrm{Total}}$ and
$\sigma^{a}_{\mathrm{Total}}$ for later use. The parameter $a$ takes on two
values: $a=0$ reduces (97) and (98) to the expressions of Anderson et al. [2],
where we set $\epsilon^{2}_{E}=\lambda_{u}$ and
$\epsilon^{2}_{S}=-\lambda_{s}$, while for $a=1$ they reduce to those used by
Guo and Lin [17] provided the choice of an isotropic surface energy is
assumed, viz., $\Gamma(\nabla c)=|\nabla c|$. Thus, as is clear from (98),
(97), and (99) that the dimensions of $\lambda_{f}$, $\lambda_{s}$, and
$\lambda_{u}$ are either specific or volumetric depending on the case. As
usual, we have the thermodynamic relation

$f=u-Ts\,,$ (100)

which allows us to assume the relationship between the coefficients

$\lambda_{f}(T)=\lambda_{u}-T\lambda_{s}\quad\text{ and
}\quad\frac{\text{d}\lambda_{f}(T)}{\text{d}T}=-\lambda_{s}\,.$ (101)

To summarize, our expressions (97) and (98) generalize the model studied by
Guo and Lin [17] by including $\Gamma$, which accounts for anisotropic surface
energy effects, while our expressions generalize the model of Anderson et al.
[2] by including the factors of $\rho$ in the integrands making all quantities
in the integrands specific quantities multiplied by the density, giving rise
to more general sources of energy.

Because $\Gamma$ is a homogeneous function of degree unity,

$\Gamma(\lambda\mathbf{p})=\lambda\Gamma(\mathbf{p})\text{ for all
}\lambda>0\,.$ (102)

Differentiating (102) with respect to $\lambda$ and then setting $\lambda=1$
yields the fundamental relation

$\Gamma(\mathbf{p})=\mathbf{p}\cdot\bm{\xi}\coloneqq
p_{j}\frac{\partial\Gamma(\mathbf{p})}{\partial p_{j}}\,.$ (103)

Then, differentiating (103) gives a second well-known relation,

$\frac{\partial\Gamma}{\partial p_{i}}=\frac{\partial}{\partial
p_{i}}(\bm{\xi}\cdot\mathbf{p})=\xi_{i}+\frac{\partial^{2}\Gamma}{\partial
p_{i}\partial p_{j}}p_{j}=\xi_{i}\,,$ (104)

where evidently $p_{j}$ must be a null eigenvector of the matrix
${\partial^{2}\Gamma}/{\partial p_{i}\partial p_{j}}$. Henceforth we will
assume the argument of $\Gamma$ to be $\nabla c$. For the case of isotropic
surface energy, where $\Gamma(\nabla c)=|\nabla c|$, the associated
homogeneous function of degree zero is given by

${\bm{\xi}}=\nabla c/|\nabla c|\,.$ (105)

From (99) we can obtain a generalized chemical potential

$\displaystyle\mu^{a}_{\Gamma}$
$\displaystyle\coloneqq\frac{\delta\mathfrak{F}^{a}}{\delta\tilde{c}}=\rho\frac{\partial
u}{\partial\tilde{c}}-\nabla\cdot(\lambda_{f}\rho^{a}\Gamma\nabla c)$
$\displaystyle=\mu-\frac{1}{\rho}\nabla\cdot(\lambda_{f}\rho^{a}\Gamma\bm{\xi})\,,$
(106)

where recall $\tilde{c}=\rho c$. For the case of isotropic surface energy,
this becomes

$\mu^{a}_{|\nabla c|}=\mu-\frac{1}{\rho}\nabla\cdot(\lambda_{f}\rho^{a}\nabla
c)\,.$ (107)

Upon setting $a=1$ [the case of 17], this reduces to

$\mu^{1}_{|\nabla c|}=\mu-\frac{1}{\rho}\nabla\cdot(\lambda_{f}\rho\,\nabla
c)\,,$ (108)

an expression that differs from that in Guo and Lin [17] unless
$\lambda_{f}\rho$ is constant. If this is the case and we choose a classical
$\mu=c^{3}-c$, corresponding to the quartic Laudau potential, we obtain

$\mu_{CH}=c^{3}-c-{\lambda_{f}}\nabla^{2}c\,,$ (109)

the chemical potential of Cahn and Hilliard who indeed make these assumptions
[cf. page 267 of 4]. For $\alpha=0$ [the case of 2], we have

$\mu^{0}_{\Gamma}=\mu-\frac{1}{\rho}\nabla\cdot(\lambda_{f}\,\Gamma\bm{\xi})\,,$
(110)

which allows for the weighted mean curvature effects of anisotropy.

Maintaining the same set of dynamical variables as in (51) (or an equivalent
set) and making the choices of energy and entropy functionals of (97) and
(98), we have completed the first two steps of the algorithm.

###  3.2 Noncanonical Poisson bracket of the Cahn-Hilliard-Euler system

To complete the next step of the algorithm, the third step, we need to
manufacture a bracket that has (98) as a Casimir invariant. We do this by
starting from the bracket of (15) in terms of the original variables
$\Psi=\\{\mathbf{m},\rho,\tilde{c},\sigma\\}$ and then transforming it to a
new set of dynamical variables
$\hat{\Psi}^{a}\\!\coloneqq\\{\hat{\mathbf{m}},\hat{\rho},\hat{\tilde{c}},\hat{\sigma}^{a}\\}$,
giving the same bracket in terms of new coordinates. We have included the
superscript $a$ because in effect we have two sets of coordinates,
corresponding to the desired entropies of (98) for $a=0$ and $a=1$. To
distinguish the old from the new, we write the bracket in the transformed
variables as $\\{\hat{F},\hat{G}\\}^{a}$. Because of coordinate invariance,
$\\{C,F\\}=\\{C^{a},\hat{F}\\}^{a}=0$, where $F[\Psi]=\hat{F}[\hat{\Psi}^{a}]$
is any functional written in one or the other coordinates. The Casimir
$S=\int_{\Omega}\sigma$ in our original coordinates is transformed into a
different form in the new coordinates. Specifically, we change the variables
as follows:

$\displaystyle{\mathbf{m}}$
$\displaystyle=\hat{\mathbf{m}}\,,\quad{\rho}=\hat{\rho}\,,\quad{\tilde{c}}=\hat{\tilde{c}}\,,$
$\displaystyle\sigma$
$\displaystyle=\hat{\sigma}^{a}+\frac{\hat{\rho}^{a}}{2}\lambda_{s}\Gamma^{2}(\nabla\hat{c})\,,$
(111)

where $\hat{\tilde{c}}=\hat{\rho}\hat{c}$. Consequently, the entropy $S$ in
the old coordinates written in terms of the new coordinates will, by design,
become the following Casimir for the Poisson bracket in the new coordinates:

$\hat{S}^{a}=\int_{\Omega}\hat{\sigma}^{a}+\frac{\hat{\rho}^{a}}{2}\lambda_{s}\Gamma^{2}(\nabla\hat{c})\,.$
(112)

Thus we have manufactured a bracket with the entropy expression of (98) as a
Casimir.

Transformation of the Poisson bracket (15) requires use of the functional
chain rule. For convenience we use

$\hat{\sigma}^{a}=\sigma-\frac{\rho^{a}}{2}\lambda_{s}\Gamma^{2}(\nabla c)$
(113)

and consider the variation of any functional of the new variables. Thus we use
$\delta\rho=\delta\hat{\rho}$, $\delta\mathbf{m}=\delta\hat{\mathbf{m}}$,
$\delta\tilde{c}=\delta\hat{\tilde{c}}$, and for the entropy variable

$\displaystyle\delta\hat{\sigma}^{a}$
$\displaystyle=\delta\sigma-\frac{1}{2}a\rho^{a-1}\lambda_{s}\Gamma^{2}(\nabla
c)\delta\rho-\rho^{a}\lambda_{s}\Gamma{\bm{\xi}}\cdot\nabla\delta\left(\frac{\tilde{c}}{\rho}\right)$
$\displaystyle=\delta\sigma-\frac{1}{2}a\rho^{a-1}\lambda_{s}\Gamma^{2}(\nabla
c)\delta\rho-\rho^{a}\lambda_{s}\Gamma{\bm{\xi}}\cdot\nabla\left(\frac{\delta\tilde{c}}{\rho}\right)$
$\displaystyle\hskip
56.9055pt+\rho^{a}\lambda_{s}\Gamma{\bm{\xi}}\cdot\nabla\left(\frac{\tilde{c}}{\rho^{2}}\delta\rho\right)\,,$
(114)

where use has been made of (104). Now let $F$ be an arbitrary functional of
the original variables and $\hat{F}$ the same functional in terms of the new
variables. Thus,

$\displaystyle\int_{\Omega}{\hat{F}}_{\hat{\mathbf{m}}}\cdot\,\delta\hat{\mathbf{m}}$
$\displaystyle+\hat{F}_{\hat{\rho}}\,{\delta\hat{\rho}}+\hat{F}_{\hat{\sigma}^{a}}\,\delta{\hat{\sigma}^{a}}+\hat{F}_{\hat{\tilde{c}}}\,\delta\hat{\tilde{c}}$
(115)
$\displaystyle=\int_{\Omega}F_{\mathbf{m}}\cdot\delta\mathbf{m}+F_{\rho}\,\delta{\rho}+F_{\sigma}\,\delta{\sigma}+F_{\tilde{c}}\,\delta{\tilde{c}}\,.$

Note, no sum over $a$ is to be done. By identification of terms we obtain

$\displaystyle F_{\mathbf{m}}$
$\displaystyle=\hat{F}_{\hat{\mathbf{m}}}\,,\quad
F_{\sigma}=\hat{F}_{\hat{\sigma}_{a}}\,,$ $\displaystyle F_{\rho}$
$\displaystyle=\hat{F}_{\hat{\rho}}-\frac{a}{2}\hat{\rho}^{a-1}\lambda_{s}\Gamma^{2}\,\hat{F}_{\hat{\sigma}^{a}}-\frac{\hat{\tilde{c}}}{\hat{\rho}^{2}}\,\nabla\cdot\left(\hat{\rho}^{a}\lambda_{s}\Gamma{\bm{\xi}}\,\hat{F}_{\hat{\sigma}^{a}}\right)\,,$
$\displaystyle F_{\tilde{c}}$
$\displaystyle=\hat{F}_{\hat{\tilde{c}}}+\frac{1}{\hat{\rho}}\nabla\cdot\left(\hat{\rho}^{a}\lambda_{s}\Gamma{\bm{\xi}}\hat{F}_{\hat{\sigma}^{a}}\right)\,.$
(116)

The transformed Poisson bracket is obtained by inserting the expressions of
(116) into (15), writing it entirely in terms of the hat variables. Upon doing
this and then dropping the hats, we get for any functionals ${F}$ and ${G}$
the following Poisson bracket:

$\displaystyle\\{{F},G\\}^{a}$
$\displaystyle=-\int_{\Omega}{\mathbf{m}}\cdot\Big{[}F_{\mathbf{m}}\cdot\nabla
G_{\mathbf{m}}-G_{\mathbf{m}}\cdot\nabla F_{\mathbf{m}}\Big{]}$
$\displaystyle+{\rho}\,\Big{[}F_{\mathbf{m}}\cdot\nabla\big{(}G_{\rho}-{a}\rho^{a-1}\lambda_{s}\Gamma^{2}G_{\sigma^{a}}/2$
$\displaystyle\hskip
56.9055pt-{\tilde{c}}\,\nabla\cdot(\rho^{a}\lambda_{s}\Gamma{\bm{\xi}}G_{\sigma^{a}})/\rho^{2}\big{)}$
$\displaystyle\hskip 14.22636pt-
G_{\mathbf{m}}\cdot\nabla\big{(}F_{\rho}-{a}\rho^{a-1}\lambda_{s}\Gamma^{2}F_{\sigma^{a}}/2$
(117) $\displaystyle\hskip
56.9055pt-{\tilde{c}}\,\nabla\cdot\left(\rho^{a}\lambda_{s}\Gamma{\bm{\xi}}F_{\sigma^{a}}\right)/\rho^{2}\big{)}\Big{]}$
$\displaystyle+\left({\sigma}^{a}+{{\rho}^{a}}\lambda_{s}\Gamma^{2}/2\right)\Big{[}F_{\mathbf{m}}\cdot\nabla
G_{\sigma^{a}}-G_{\mathbf{m}}\cdot\nabla F_{\sigma^{a}}\Big{]}$
$\displaystyle+\tilde{c}\,\Big{[}F_{\mathbf{m}}\cdot\nabla\left(G_{\tilde{c}}+\nabla\cdot\left(\rho^{a}\lambda_{s}\Gamma{\bm{\xi}}G_{\sigma^{a}}\right)/\rho\right)$
$\displaystyle\hskip 14.22636pt-
G_{\mathbf{m}}\cdot\nabla\left(F_{\tilde{c}}+\nabla\cdot\left(\rho^{a}\lambda_{s}\Gamma{\bm{\xi}}F_{\sigma^{a}}\right)/\rho\right)\Big{]}\,.$

This bracket is clearly bilinear and skew-symmetric. Because it was derived
from the bracket 15 by a change of variables, satisfaction of the Jacobi
identity is assured. We note, as before, strong boundary conditions are
assumed such that all integrations by parts produce vanishing boundary terms.

Thus we have completed the third part of our algorithm, the construction of a
Poisson bracket that has entropy functional of (98) in the set of its Casimir
invariants. Recall the integrand of the entropy is given by

$\sigma_{\mathrm{Total}}^{a}\coloneqq\sigma^{a}+\frac{\rho^{a}}{2}\lambda_{s}\Gamma^{2}(\nabla
c)\,;$ (118)

so we find

$\displaystyle\frac{\delta S^{a}}{\delta\sigma^{a}}$
$\displaystyle=1\,,\quad\frac{\delta
S^{a}}{\delta\tilde{c}}=-\frac{1}{\rho}\nabla\cdot\left(\rho^{a}\lambda_{s}\Gamma{\bm{\xi}}\right)$
$\displaystyle\frac{\delta S^{a}}{\delta\rho}$
$\displaystyle=\frac{a}{2}\rho^{a-1}\lambda_{s}\Gamma^{2}+\frac{\tilde{c}}{\rho^{2}}\nabla\cdot\left(\rho^{a}\lambda_{s}\Gamma{\bm{\xi}}\right)\,.$
(119)

Using (119) one can easily check that $\\{F,S^{a}\\}=0$ for all $F$, which by
construction had to be the case. Now we are free to choose any Hamiltonian we
desire in (117) to obtain the evolution of any observable $o$ as follows:

$\partial_{t}o=\\{o,H^{a}\\}^{a}\,,$ (120)

In §§ 3.1 we proposed the Hamiltonian functional of (97), which we rewrite as
follows in order to make all arguments clear:

$\displaystyle H^{a}[\rho,\mathbf{m},\sigma,\tilde{c}]$
$\displaystyle=\int_{\Omega}\frac{|\mathbf{m}|^{2}}{2\rho}+\rho
u\left(\rho,\frac{\sigma^{a}}{\rho},\frac{\tilde{c}}{\rho}\right)$
$\displaystyle\hskip
45.52458pt+\frac{\rho^{a}}{2}\lambda_{u}\Gamma^{2}\left(\nabla\left(\frac{\tilde{c}}{\rho}\right)\right)\,.$
(121)

Using the functional derivatives of this Hamiltonian,

$\displaystyle H^{a}_{\mathbf{m}}$ $\displaystyle=\mathbf{v},\quad
H_{\sigma^{a}}=T\,,$ $\displaystyle H^{a}_{\rho}$
$\displaystyle=-|\mathbf{v}|^{2}/2+u+{p}/{\rho}-sT-c\mu$ $\displaystyle\hskip
28.45274pt+a\rho^{a-1}\lambda_{u}\Gamma^{2}/2+{c}\,\nabla\cdot(\rho^{a}\lambda_{u}\Gamma{\bm{\xi}})/\rho\,,$
$\displaystyle H^{a}_{\tilde{c}}$
$\displaystyle=\mu-\nabla\cdot(\rho^{a}\lambda_{u}\Gamma{\bm{\xi}})/\rho\,,$
(122)

in the bracket (117) gives the equations of motion in the form of (120). At
this point we could write this out and display a general system of equations
that includes both cases, but we choose to consider them separately because
the general system is unwieldy and not particularly perspicuous.

Let us first consider the simplified version of our derived Poisson bracket
for the case $a=1$, which is as follows:

$\displaystyle\\{F,G\\}^{1}$
$\displaystyle=-\int_{\Omega}\mathbf{m}\,\cdot\big{[}F_{\mathbf{m}}\cdot\nabla
G_{\mathbf{m}}-G_{\mathbf{m}}\cdot\nabla F_{\mathbf{m}}\big{]}$
$\displaystyle+\rho\big{[}F_{\mathbf{m}}\cdot\nabla
G_{\rho}-G_{\mathbf{m}}\cdot\nabla F_{\rho}\big{]}$
$\displaystyle-\lambda_{s}\big{[}F_{\mathbf{m}}\cdot\nabla\cdot\left(\rho
G_{\sigma^{1}}\Gamma{\bm{\xi}}\otimes\nabla c\right)$ $\displaystyle\hskip
56.9055pt-G_{\mathbf{m}}\cdot\nabla\cdot\left(\rho
F_{\sigma^{1}}\Gamma{\bm{\xi}}\otimes\nabla c\right)\big{]}$
$\displaystyle+\sigma^{1}\left[F_{\mathbf{m}}\cdot\nabla
G_{\sigma^{1}}-G_{\mathbf{m}}\cdot\nabla F_{\sigma^{1}}\right]$
$\displaystyle+\tilde{c}\left[F_{\mathbf{m}}\cdot\nabla
G_{\tilde{c}}-G_{\mathbf{m}}\cdot\nabla F_{\tilde{c}}\right]\,,$ (123)

where $\otimes$ denotes tensor product of two vectors and consistent with our
convention we have

$\nabla\cdot(\mathbf{u}\otimes\mathbf{v})=(\nabla\cdot\mathbf{u})\mathbf{v}+\mathbf{u}\cdot\nabla\mathbf{v}\,.$

Using (122), the bracket form of (120) gives the ideal diffuse two-phase flow
system

$\displaystyle\partial_{t}\mathbf{v}$
$\displaystyle=\\{\mathbf{v},H^{1}\\}^{1}$
$\displaystyle=-\mathbf{v}\cdot\nabla\mathbf{v}-\frac{1}{\rho}\nabla\cdot(p\bar{I}+\rho\lambda_{f}\Gamma{\bm{\xi}}\otimes\nabla
c)\,,$ (124) $\displaystyle\partial_{t}\rho$
$\displaystyle=\\{\rho,H^{1}\\}^{1}=-\mathbf{v}\cdot\nabla\rho-\rho\,\nabla\cdot\mathbf{v}\,,$
(125) $\displaystyle\partial_{t}\tilde{c}$
$\displaystyle=\\{\tilde{c},H^{1}\\}^{1}=-\mathbf{v}\cdot\nabla\tilde{c}-\tilde{c}\,\nabla\cdot\mathbf{v}\,,$
(126) $\displaystyle\partial_{t}\sigma^{1}_{\mathrm{Total}}$
$\displaystyle=\\{{\sigma}^{1}_{\mathrm{Total}},H^{1}\\}^{1}$
$\displaystyle=-\mathbf{v}\cdot\nabla{\sigma}^{1}_{\mathrm{Total}}-{\sigma}^{1}_{\mathrm{Total}}\,\nabla\cdot\mathbf{v}\,,$
(127)

where $\bar{I}$ is the unit tensor. Observe in (127) we have chosen the
observable $\sigma^{1}_{\mathrm{Total}}$ instead of $\sigma^{1}$, in order to
demonstrate its conservation.

Similarly, for the case where $a=0$, the Poisson bracket has the following
form:

$\displaystyle\\{F,G\\}^{0}$
$\displaystyle=-\int_{\Omega}\mathbf{m}\cdot\,\big{[}F_{\mathbf{m}}\cdot\nabla
G_{\mathbf{m}}-G_{\mathbf{m}}\cdot\nabla F_{\mathbf{m}}\big{]}$
$\displaystyle+\rho\big{[}F_{\mathbf{m}}\cdot\nabla
G_{\rho}-G_{\mathbf{m}}\cdot\nabla F_{\rho}\big{]}$
$\displaystyle-\lambda_{s}\big{[}F_{\mathbf{m}}\cdot\nabla\cdot\left(G_{\sigma^{0}}\Gamma{\bm{\xi}}\otimes\nabla
c\right)$ $\displaystyle\hskip 68.28644pt-
G_{\mathbf{m}}\cdot\nabla\cdot\left(F_{\sigma^{0}}\Gamma{\bm{\xi}}\otimes\nabla
c\right)\big{]}$
$\displaystyle+{\lambda_{s}}\left[F_{\mathbf{m}}\cdot\nabla\left(\Gamma^{2}G_{\sigma^{0}}\right)-G_{\mathbf{m}}\cdot\nabla\left(\Gamma^{2}F_{\sigma^{0}}\right)\right]/2$
$\displaystyle+\sigma^{0}\left[F_{\mathbf{m}}\cdot\nabla
G_{\sigma^{0}}-G_{\mathbf{m}}\cdot\nabla F_{\sigma^{0}}\right]$
$\displaystyle+\tilde{c}\left[F_{\mathbf{m}}\cdot\nabla
G_{\tilde{c}}-G_{\mathbf{m}}\cdot\nabla F_{\tilde{c}}\right]\,.$ (128)

Same as above, using (122), the ideal diffuse two-phase flow system is
produced

$\displaystyle\partial_{t}\mathbf{v}$
$\displaystyle=\\{\mathbf{v},H^{0}\\}^{0}$ (129)
$\displaystyle=-\mathbf{v}\cdot\nabla\mathbf{v}-\frac{1}{\rho}\nabla\cdot\Big{[}\left(p-\lambda_{f}\Gamma^{2}/2\right)\bar{I}$
$\displaystyle\hskip 92.47145pt+\lambda_{f}\Gamma{\bm{\xi}}\otimes\nabla
c\Big{]}\,,$ $\displaystyle\partial_{t}\rho$
$\displaystyle=\\{\rho,H^{0}\\}^{0}=-\mathbf{v}\cdot\nabla\rho-\rho\,\nabla\cdot\mathbf{v}\,,$
(130) $\displaystyle\partial_{t}\tilde{c}$
$\displaystyle=\\{\tilde{c},H^{0}\\}^{0}=-\mathbf{v}\cdot\nabla\tilde{c}-\tilde{c}\,\nabla\cdot\mathbf{v}\,,$
(131) $\displaystyle\partial_{t}\sigma^{0}_{\mathrm{Total}}$
$\displaystyle=\\{{\sigma}^{0}_{\mathrm{Total}},H^{0}\\}^{0}$
$\displaystyle=-\mathbf{v}\cdot\nabla{\sigma}_{\mathrm{Total}}^{0}-{\sigma}_{\mathrm{Total}}^{0}\,\nabla\cdot\mathbf{v}\,.$
(132)

Let us now comment on these two Hamiltonian systems. By construction both the
$a=1$ and $a=0$ systems conserve their Hamiltonians and entropies, as given by
(98) and (97) with $a=1$ and $a=0$ , respectively. Both systems have momentum
equations containing a term describing anisotropic surface energy (capillary)
effects. The $a=0$ system of (129)–(132) is identical to the ideal limit of
that given in the work of Anderson et al. [2]. Upon choosing $\Gamma(\nabla
c)=|\nabla c|$, the $a=1$ system of (124)–(127) should correspond to the ideal
limit of that of Guo and Lin [17], but it does not. In fact the system of Guo
and Lin [17] in this limit does not conserve energy. Moreover, the capillary
effect in their momentum equation (equation (3.40)), which should be replaced
by (124) with $\Gamma(\nabla c)=|\nabla c|$, vanishes in the one-dimensional
limit. Since such surface effects are determined by mean or weighted mean
curvature [34], it is clear that this is physically untenable. Fortunately,
our method provides a simple fix to their equations, while showing how to
generalize them to include anisotropic surface effects.

An alternative but equivalent Hamiltonian formulation of the above systems
exists, in fact, one that has a standard entropy functional of the form of
(30). Given that the bracket of (117) was obtained via a transformation of the
bracket of (15), we can transform it back from one that has (98) as a Casimir
to the original that has (5) as a Casimir. However, to generate equivalent
equations of motion, we would have to transform the Hamiltonian of (97) into a
more complicated form. Tracing back through our transformations, we would
replace the coordinate $\sigma^{a}$ in the Hamiltonian by
$\sigma-{\rho^{a}}\lambda_{s}\Gamma^{2}/2$, which means the internal energy
becomes

$u\mapsto
u(\rho,(\sigma-{\rho^{a}}\lambda_{s}\Gamma^{2}/2)/{\rho},\tilde{c}/\rho)\,,$
(133)

while otherwise the Hamiltonian remains the same. Just as with finite-
dimensional Hamiltonian systems, one can change coordinates and arrive at
equivalent systems with different Poisson brackets and Hamiltonians, and in
the noncanonical case different expressions for the Casimir invariants. Often
one has the options of a simple bracket and complicated Hamiltonian or vice
verse.

###  3.3 Metriplectic 4-bracket for the Cahn-Hilliard-Navier-Stokes system

Now we turn to the 4th and final step of our algorithm. Just as in §§ 2.3.3 we
build a metriplectic 4-bracket using the K-N construction. We suppose our
multi-component field variable is

$\psi(\mathbf{x},t)=(\mathbf{m}(\mathbf{x},t),\rho(\mathbf{x},t),\tilde{c}(\mathbf{x},t),\sigma^{a}(\mathbf{x},t))$
(134)

and in the K-N construction we use a more general expression for $\Sigma$ akin
to that mentioned in (56), viz.

$\displaystyle M(dF,dG)$
$\displaystyle=F_{{{\sigma}^{a}}}G_{{{\sigma}^{a}}}\,,$ (135)
$\displaystyle{\Sigma}(dF,dG)$ $\displaystyle=\nabla
F_{{\mathbf{m}}}:\bar{\bar{\Lambda}}_{1}:\nabla G_{{\mathbf{m}}}+\nabla
F_{{{\sigma}^{a}}}\cdot\bar{\Lambda}_{2}\cdot\nabla G_{{{\sigma}^{a}}}$
$\displaystyle\hskip
42.67912pt+\nabla\mathcal{L}^{a}_{\tilde{c}}(F)\cdot\bar{\Lambda}_{3}\cdot\mathcal{L}^{a}_{\tilde{c}}(F)\,,$
(136)

where $a$, of course, is not to be summed over and the pseudodifferential
operator $\mathcal{L}^{a}_{\tilde{c}}$ has the following form:

$\mathcal{L}^{a}_{\tilde{c}}(F):=\nabla\big{(}F_{\tilde{c}}+\nabla\cdot\left(\rho^{a}\lambda_{s}\Gamma{\bm{\xi}}F_{\sigma^{a}}\right)/\rho\big{)}\,.$
(137)

Here the tensors $\bar{\bar{\Lambda}}_{1}$ , $\bar{\Lambda}_{2}$ and
$\bar{\Lambda}_{3}$ are defined by (59). Then, the 4-bracket reads

$\displaystyle(F,K;G,N)^{a}$ $\displaystyle=$ (138)
$\displaystyle\int_{\Omega}\\!\frac{1}{T}\Big{[}\left[K_{{\sigma}^{a}}\nabla
F_{\mathbf{m}}-F_{{\sigma}^{a}}\nabla
K_{\mathbf{m}}\right]\\!\colon\\!\\!\bar{\bar{\Lambda}}\colon\\!\\!\left[N_{{\sigma}^{a}}\nabla
G_{\mathbf{m}}-G_{{\sigma}^{a}}\nabla N_{\mathbf{m}}\right]$
$\displaystyle+\frac{1}{T}\big{[}K_{{\sigma}^{a}}\nabla
F_{{{\sigma}^{a}}}-F_{{\sigma}^{a}}\nabla
K_{{{\sigma}^{a}}}\big{]}\cdot\bar{\kappa}\cdot\big{[}N_{{\sigma}^{a}}\nabla
G_{{{\sigma}^{a}}}-G_{{\sigma}^{a}}\nabla N_{{{\sigma}^{a}}}\big{]}$
$\displaystyle+\\!\big{[}K_{{\sigma}^{a}}\mathcal{L}^{a}_{\tilde{c}}(F)-F_{{\sigma}^{a}}\mathcal{L}^{a}_{\tilde{c}}(K)\big{]}\\!\cdot\\!\bar{D}\\!\cdot\\!\big{[}N_{{\sigma}^{a}}\mathcal{L}^{a}_{\tilde{c}}(G)-G_{{\sigma}^{a}}\mathcal{L}^{a}_{\tilde{c}}(N)\big{]}.$

Observe, with the exception of the last line this bracket is identical to that
of (60).

Upon insertion of $H^{a}$ as given by (97) and $S$ from the set of Casimirs to
be as in (98), the dynamics is given by

$\partial_{t}\psi^{\alpha}=\\{\psi^{\alpha},H^{a}\\}^{a}+(\psi^{\alpha},H^{a};S^{a},H^{a})^{a}\,.$
(139)

Using $\mathcal{L}^{a}_{\tilde{c}}({H^{a}})=\nabla\mu^{a}_{\Gamma}$,
$H^{a}_{\mathbf{m}}=\mathbf{v}$, $H^{a}_{{\sigma}_{a}}=T$,
$S^{a}_{{\sigma}_{a}}=1$ and $\mathcal{L}^{a}_{\tilde{c}}({S^{a}})=0$, the
following diffuse-interface CHNS system for $a=1$ is produced:

$\displaystyle\partial_{t}\mathbf{v}$
$\displaystyle=\\{\mathbf{v},H^{1}\\}^{1}+(\mathbf{v},H^{1};S^{1},H^{1})^{1}$
$\displaystyle=-\mathbf{v}\cdot\nabla\mathbf{v}-\frac{1}{\rho}\nabla\cdot\big{[}p\mathbf{I}+\lambda_{f}\rho\Gamma{\bm{\xi}}\otimes\nabla
c\big{]}$ $\displaystyle\hskip
85.35826pt+\frac{1}{\rho}\nabla\cdot(\bar{\bar{\Lambda}}:\nabla\mathbf{v})\,,$
(140) $\displaystyle\partial_{t}\rho$
$\displaystyle=\\{\rho,H^{1}\\}^{1}+(\rho,H^{1};S^{1},H^{1})^{1}$
$\displaystyle=-\mathbf{v}\cdot\nabla\rho-\rho\,\nabla\cdot\mathbf{v}\,,$
(141) $\displaystyle\partial_{t}\tilde{c}$
$\displaystyle=\\{\tilde{c},H^{1}\\}^{1}+(\tilde{c},H^{1};S^{1},H^{1})^{1}$
$\displaystyle=-\mathbf{v}\cdot\nabla\tilde{c}-\tilde{c}\,\nabla\cdot\mathbf{v}+\nabla\cdot(\bar{D}\cdot\nabla\mu^{1}_{\Gamma})\,,$
(142) $\displaystyle\partial_{t}\sigma^{1}_{\mathrm{Total}}$
$\displaystyle=\\{\sigma^{1}_{\mathrm{Total}},H^{1}\\}^{1}+(\sigma^{1}_{\mathrm{Total}},H^{1};S^{1},H^{1})^{1}$
$\displaystyle=-\mathbf{v}\cdot\nabla\sigma^{1}_{\mathrm{Total}}-\sigma^{1}_{\mathrm{Total}}\,\nabla\cdot\mathbf{v}$
$\displaystyle\quad+\nabla\cdot\left(\frac{\bar{\kappa}}{T}\cdot\nabla
T\right)+\frac{1}{T^{2}}\nabla T\cdot\bar{\kappa}\cdot\nabla T$ (143)
$\displaystyle\hskip
28.45274pt+\frac{1}{T}\nabla\mathbf{v}:\bar{\bar{\Lambda}}:\nabla\mathbf{v}+\frac{1}{T}\nabla\mu^{1}_{\Gamma}\cdot\bar{D}\cdot\nabla\mu^{1}_{\Gamma}\,.$

Similarly, for $a=0$ we obtain

$\displaystyle\partial_{t}\mathbf{v}$
$\displaystyle=\\{\mathbf{v},H^{0}\\}^{0}+(\mathbf{v},H^{0};S^{0},H^{0})^{0}$
$\displaystyle=-\mathbf{v}\cdot\nabla\mathbf{v}-\frac{1}{\rho}\nabla\cdot\Big{[}\left(p-\lambda_{f}\Gamma^{2}/2\right)\mathbf{I}$
$\displaystyle\hskip 28.45274pt+\lambda_{f}\Gamma{\bm{\xi}}\otimes\nabla
c\Big{]}+\frac{1}{\rho}\nabla\cdot(\bar{\bar{\Lambda}}:\nabla\mathbf{v})\,,$
(144) $\displaystyle\partial_{t}\rho$
$\displaystyle=\\{\rho,H^{0}\\}^{0}+(\rho,H^{0};S^{0},H^{0})^{0}$
$\displaystyle=-\mathbf{v}\cdot\nabla\rho-\rho\,\nabla\cdot\mathbf{v}\,,$
(145) $\displaystyle\partial_{t}\tilde{c}$
$\displaystyle=\\{\tilde{c},H^{0}\\}^{0}+(\tilde{c},H^{0};S^{0},H^{0})^{0}$
$\displaystyle=-\mathbf{v}\cdot\nabla\tilde{c}-\tilde{c}\,\nabla\cdot\mathbf{v}+\nabla\cdot(\bar{D}\cdot\nabla\mu^{0}_{\Gamma})\,,$
(146) $\displaystyle\partial_{t}\sigma^{0}_{\mathrm{Total}}$
$\displaystyle=\\{\sigma^{0}_{\mathrm{Total}},H^{0}\\}^{0}+(\sigma^{0}_{\mathrm{Total}},H^{0};S^{0},H^{0})^{0}$
$\displaystyle=-\mathbf{v}\cdot\nabla\sigma^{0}_{\mathrm{Total}}-\sigma^{0}_{\mathrm{Total}}\,\nabla\cdot\mathbf{v}$
(147) $\displaystyle\quad+\nabla\cdot\left(\frac{\bar{\kappa}}{T}\cdot\nabla
T\right)+\frac{1}{T^{2}}\nabla T\cdot\bar{\kappa}\cdot\nabla T$
$\displaystyle\hskip
28.45274pt+\frac{1}{T}\nabla\mathbf{v}:\bar{\bar{\Lambda}}:\nabla\mathbf{v}+\frac{1}{T}\nabla\mu^{0}_{\Gamma}\cdot\bar{D}\cdot\nabla\mu^{0}_{\Gamma}\,.$

Thus we have extracted from our general system with arbitrary $a$, two
thermodynamically consistent CHNS systems. By construction both the $a=1$ and
$a=0$ systems must conserve energy and both must produce entropy, which we
find is governed by the following:

$\displaystyle\dot{S^{a}}$ $\displaystyle=(S^{a},H^{a};S^{a},H^{a})^{a}$
$\displaystyle=\int_{\Omega}\frac{1}{T}\bigg{[}\nabla\mathbf{v}:\bar{\bar{\Lambda}}:\nabla\mathbf{v}+\frac{1}{T}\nabla
T\cdot\bar{\kappa}\cdot\nabla T$ $\displaystyle\hskip
88.2037pt+\nabla\mu^{a}_{\Gamma}\cdot\bar{D}\cdot\nabla\mu^{a}_{\Gamma}\bigg{]}\geq
0\,.$ (148)

###  3.4 Metriplectic 2-bracket for the Cahn-Hilliard-Navier-Stokes system

Proceeding as in the previous section, using
$\mathcal{L}^{a}_{\tilde{c}}({H}^{a})=\nabla(\mu-\frac{1}{\rho}\nabla\cdot(\lambda_{f}\rho^{a}\Gamma\bm{\xi}))=\nabla\mu_{\Gamma}^{a}$,
$H^{a}_{\mathbf{m}}=\mathbf{v}$, $H^{a}_{{\sigma}^{a}}=T$, the metriplectic
2-bracket emerges directly from the 4-bracket as follows:

$\displaystyle(F,{G})^{a}_{{H^{a}}}$
$\displaystyle=(F,{H^{a}};{G},{H^{a}})^{a}$ (149)
$\displaystyle=\int_{\Omega}\frac{1}{T}\Big{[}\left[T\nabla
F_{\mathbf{m}}-F_{{\sigma}^{a}}\nabla\mathbf{v}\right]:\bar{\bar{\Lambda}}:\left[T\nabla
G_{\mathbf{m}}-G_{{\sigma}^{a}}\nabla\mathbf{v}\right]$ $\displaystyle+\
\frac{1}{T}\,\big{[}T\nabla F_{{{\sigma}^{a}}}-F_{{\sigma}^{a}}\nabla
T\big{]}\cdot\bar{\kappa}\cdot\big{[}T\nabla
G_{{{\sigma}^{a}}}-G_{{\sigma}^{a}}\nabla T\big{]}$
$\displaystyle+\big{[}T\mathcal{L}^{a}_{\tilde{c}}(F)-F_{\sigma^{a}}\nabla\mu_{\Gamma}^{a}\big{]}\\!\cdot\bar{D}\cdot\\!\big{[}T\mathcal{L}^{a}_{\tilde{c}}(G)-G_{\sigma^{a}}\nabla\mu_{\Gamma}^{a}\big{]}\Big{]}\,,$

which is the analog of (60) for the GNS system in §§ 2.3.3.

Now we could proceed as in §§ 2.4.1 and obtain the analog of $L_{\alpha\beta}$
by transforming to the coordinates
$\xi^{a}\\!\\!\coloneqq(\mathbf{m},\rho,\tilde{c},e^{a}_{\mathrm{Total}})$,
where $e^{a}_{\mathrm{Total}}$ is the total energy density defined in (97).
This would lead to a metriplectic 2-bracket analogous to (76) and the
concomitant flux and affinity relations analogous to (87) would emerge.
Because this is hardly more enlightening than (149), we do not record this
result here.

##  4 Summary and Conclusions

In this paper we have described the metriplectic 4-bracket formalism and how
it can algorithmically be used in the context of multiphase fluids to
construct thermodynamically consistent models, ones that conserve energy and
produce entropy. In particular, we have used it in § 2 to obtain the GE
Hamiltonian system, which adds a concentration variable to conventional
Eulerian fluid mechanics, and to the GNS system, a generalization of the
Navier-Stokes system that is thermodynamically consistent with the collection
of thermodynamic fluxes. Then, in § 3 we used the algorithm to obtain a class
of Hamiltonian fluid systems that allow for anisotropic surface effects,
followed by the construction of a general class of CHNS systems that couple
Cahn-Hilliard physics with that of Navier-Stokes dynamics, in a
thermodynamically consistent way. The systems we obtain generalizes previous
work by including anisotropic effects in the surface tension and all
phenomenological parameters.

A cornerstone of Hamiltonian dynamics is its geometric invariance under
coordinate changes. Because the minimal metriplectic properties are algebraic
and geometric, they too are invariant under coordinate changes. Thus, we can
write our CHNS class of dissipative systems with a standard entropy functional
of the form of (30), but with a more complicated Hamiltonian using (133).

From the examples presented, it is clear that the 4-bracket formalism can be
applied to obtain a wide variety of dynamical systems in various fields. In
fact it was recently applied to obtain generalized collision operators in
kinetic theory [32] and a thermodynamically consistent model for radiation
hydrodynamics [37].

The metriplectic 4-bracket formalism also provides an avenue for designing
structure preserving numerical algorithms [see e.g. 27]. Any discretization
that preserves the symmetries of the 4-bracket, which is not a difficult task,
will be thermodynamically consistent on the semi-discrete level, i.e. produce
a set of ordinary differential equations that conserve energy and produce
entropy.

##  Appendix A Semidefinite curvature

We define the binary operations $\langle\cdot,\cdot\rangle_{\Sigma}$ and
$\langle\cdot,\cdot\rangle_{M}$ that satisfy all of the axioms of an inner
product space, except the non-degeneracy condition

$\displaystyle\langle F,G\rangle_{\Sigma}$
$\displaystyle:=\int\\!\\!d^{N}\\!z\\!\\!\int\\!\\!d^{N}\\!z^{\prime}\,\Sigma^{\alpha\beta}(z,z^{\prime})\frac{\delta
F}{\delta\chi^{\alpha}(z)}\frac{\delta G}{\delta\chi^{\beta}(z^{\prime})}\,,$
$\displaystyle\langle F,G\rangle_{M}$
$\displaystyle:=\\!\\!\int\\!\\!d^{N}\\!z^{\prime\prime}\\!\\!\int\\!\\!d^{N}\\!z^{\prime\prime\prime}\,M^{\gamma\delta}(z^{\prime\prime},z^{\prime\prime\prime})\frac{\delta
F}{\delta\chi^{\gamma}(z^{\prime\prime})}\frac{\delta
G}{\delta\chi^{\delta}(z^{\prime\prime\prime})}\,,$

where $\Sigma$ and $M$ are positive semi-definites. We have the Cauchy-Schwarz
inequality

$\left|\langle F,G\rangle_{\Sigma}\right|\leq\sqrt{\langle
F,F\rangle_{\Sigma}}\,\,\sqrt{\langle
G,G\rangle_{\Sigma}}=\|F\|_{\Sigma}\|G\|_{\Sigma}\,.$

###### Lemma 1

A metriplectic quadravector constructed using the $K-N$ product, has non-
negative sectional curvature,

$\displaystyle K(F,G)$
$\displaystyle=\\!\\!\\!\int\\!\\!d^{N}\\!z\\!\\!\int\\!\\!d^{N}\\!z^{\prime}\,\\!\\!\\!\int\\!\\!d^{N}\\!z^{\prime\prime}\\!\\!\\!\int\\!\\!d^{N}\\!z^{\prime\prime\prime}\,\Sigma^{ij}(z,z^{\prime})\,M^{kl}(z^{\prime\prime}\\!\\!,z^{\prime\prime\prime})$
$\displaystyle\times\ \frac{\delta F}{\delta\chi^{\alpha}(z)}\frac{\delta
G}{\delta\chi^{\beta}(z^{\prime})}\,\frac{\delta
F}{\delta\chi^{\gamma}(z^{\prime\prime})}\frac{\delta
G}{\delta\chi^{\delta}(z^{\prime\prime\prime})}\ +\mathrm{\ other\ terms}\,.$

###### Proof 1

Direct calculation gives

$K(F,G)=\|F\|_{\Sigma}^{2}\|G\|_{M}^{2}-2\langle F,G\rangle_{\Sigma}\langle
F,G\rangle_{M}+\|G\|_{\Sigma}^{2}\|F\|_{M}^{2}\,.$

The following inequality

$\left(\|F\|_{\Sigma}\|G\|_{M}-\|G\|_{\Sigma}\|F\|_{M}\right)^{2}\geq 0$

implies

$\displaystyle\|F\|_{\Sigma}^{2}\|G\|_{M}^{2}+\|G\|_{\Sigma}^{2}\|F\|_{M}^{2}$
$\displaystyle\geq 2\|F\|_{M}\|F\|_{\Sigma}\|G\|_{M}\|G\|_{\Sigma}$
$\displaystyle\geq 2\left|\langle F,G\rangle_{\Sigma}\right|\left|\langle
F,G\rangle_{M}\right|$ $\displaystyle\geq 2\langle F,G\rangle_{\Sigma}\langle
F,G\rangle_{M}\,,$

where the second inequality follows from the Cauchy-Schwarz inequality.
Evidently, the last inequality implies $K(F,G)\geq 0$ for all $F$ and $G$.

###### Lemma 2

We suppose that $\Sigma$ is positive definite, defining an inner product.
Given any two $\Sigma$-arbitrary linearly independent $\delta F/\delta{\chi}$
and $\delta G/\delta{\chi}$, then the sectional curvature is strictly positive
($K(F,G)>0$).

###### Proof 2

Since $\delta F/\delta{\chi}$ and $\delta G/\delta{\chi}$ are
$\Sigma$-Linearly independent, the Cauchy-Schwarz inequality given by

$\left|\langle F,G\rangle_{\Sigma}\right|<\|F\|_{\Sigma}\|G\|_{\Sigma}\,.$

In the same way we have

$\left(\|F\|_{\Sigma}\|G\|_{M}-\|G\|_{\Sigma}\|F\|_{M}\right)^{2}\geq 0$

implies

$\displaystyle\|F\|_{\Sigma}^{2}\|G\|_{M}^{2}+\|G\|_{\Sigma}^{2}\|F\|_{M}^{2}$
$\displaystyle\geq 2\|F\|_{M}\|F\|_{\Sigma}\|G\|_{M}\|G\|_{\Sigma}$
$\displaystyle>2\left|\langle F,G\rangle_{\Sigma}\right|\left|\langle
F,G\rangle_{M}\right|$ $\displaystyle>2\langle F,G\rangle_{\Sigma}\langle
F,G\rangle_{M}\,.$

Hence, we deduce that $K(F,G)>0$ .

Finite-dimensional versions of these two lemmas were first reported in
Morrison and Updike [30].

## Acknowledgements

A.Z. acknowledges support from the Mohammed VI Polytechnic University for
supporting an internship at the University of Texas at Austin. P.J.M.
acknowledges support from the DOE Office of Fusion Energy Sciences under DE-
FG02-04ER-54742, and from a Forschungspreis from the Alexander von Humboldt
Foundation.

## References

  * Abdelhamid et al. [2015] Abdelhamid, H.M., Kawazura, Y., Yoshida, Z., 2015.  Hamiltonian formalism of extended magnetohydrodynamics.  J. Phys. A 48, 235502\.
  * Anderson et al. [2000] Anderson, D.M., McFadden, G.B., Wheeler, G.B., 2000.  A phase-field model of solidification with convection.  Physica D 135, 175–194.
  * Beris and Edwards [1994] Beris, A.N., Edwards, B.J., 1994\.  Thermodynamics of Flowing Systems with Internal Microstructure.  Number 36 in Oxford Engineering Science Series, Oxford University Press, Oxford, New York.
  * Cahn and Hilliard [1958] Cahn, J.W., Hilliard, J.E., 1958\.  Free energy of a nonuniform system. i. interfacial free energy.  J. Chem. Phys 28, 258–267.
  * Callen [1966] Callen, H.B., 1966.  Thermodynamics.  John Wiley & Sons Inc., New York.
  * Coquinot and Morrison [2020] Coquinot, B., Morrison, P.J., 2020\.  A general metriplectic framework with application to dissipative extended magnetohydrodynamics.  J. Plasma Phys. 86, 835860302\.
  * D’Avignon et al. [2016] D’Avignon, E.C., Morrison, P.J., Lingam, M., 2016.  Derivation of the Hall and extended magnetohydrodynamics brackets.  Phys. Plasmas 23, 062101\.
  * Eckart [1940a] Eckart, C., 1940a.  The thermodynamics of irreversible processes i. the simple fluid.  Phys. Rev. 58, 267–269.
  * Eckart [1940b] Eckart, C., 1940b.  The thermodynamics of irreversible processes ii. fluid mixtures.  Phys. Rev. 58, 267–269.
  * Edwards [1998] Edwards, B.J., 1998.  An analysis of single and double generator thermodynamic formalisms for the macroscopic description of complex fluids.  Journal of Non-Equilibrium Thermodynamics 23, 301–333.
  * ten Eikelder et al. [2023] ten Eikelder, M.F.P., Akkerman, I., Schillinger, D., 2023.  A unified framework for Navier-Stokes Cahn–Hilliard models with non-matching densities.  Math. Models Methods Appl. Sci. 33, 175–221.
  * Eldred and Gay-Balmaz [2020] Eldred, C., Gay-Balmaz, F., 2020\.  Single and double generator bracket formulations of multicomponent fluids with irreversible processes.  J. Phys. A 53, 395701\.
  * Fiedler [2003] Fiedler, B., 2003.  Determination of the structure algebraic curvature tensors by means of Young symmetrizers.  Séminair. Lotharingien de Combinatoire 48, B48d.
  * Grmela [1984] Grmela, M., 1984.  Particle and bracket formulations of kinetic equations.  Contemp. Math. 28, 125–132.
  * Grmela and Öttinger [1997] Grmela, M., Öttinger, H.C., 1997\.  Dynamics and thermodynamics of complex fluids. part i. development of a general formalism.  Phys. Rev. E 56, 6620–6632.
  * de Groot and Mazur [1962] de Groot, S.R., Mazur, P., 1962\.  Non-equilibrium thermodynamics.  North Holland Publishing Company, Amsterdam.
  * Guo and Lin [2015] Guo, Z., Lin, P., 2015.  A thermodynamically consistent phase-field model for two-phase flows with thermocapillary effects.  J. Fluid Mech. 766, 226–271.  doi:10.1017/jfm.2014.696.
  * Kaufman [1984] Kaufman, A.N., 1984.  Dissipative Hamiltonian systems: A unifying principle.  Phys. Lett. A 100, 419–422.
  * Kaufman and Morrison [1982] Kaufman, A.N., Morrison, P.J., 1982\.  Algebraic structure of the plasma quasilinear equations.  Phys. Lett. A 88, 405–406.
  * Kulkarni [1972] Kulkarni, R.S., 1972.  On the Bianchi identities.  Math. Ann 199, 175–204.
  * Morrison [1982] Morrison, P.J., 1982.  Poisson Brackets for Fluids and Plasmas.  AIP Conf. Proc. 88, 13–46.
  * Morrison [1984a] Morrison, P.J., 1984a.  Bracket formulation for irreversible classical fields.  Phys. Lett. A 100, 423–427.
  * Morrison [1984b] Morrison, P.J., 1984b.  Some Observations Regarding Brackets and Dissipation.  Technical Report PAM–228. University of California at Berkeley.
  * Morrison [1986] Morrison, P.J., 1986.  A paradigm for joined Hamiltonian and dissipative systems.  Physica D 18, 410–419.
  * Morrison [1998] Morrison, P.J., 1998.  Hamiltonian description of the ideal fluid.  Rev. Mod. Phys. 70, 467–521.
  * Morrison [2009] Morrison, P.J., 2009.  Thoughts on brackets and dissipation: Old and new.  J. Physics: Conf. Ser. 169, 012006\.
  * Morrison [2017] Morrison, P.J., 2017.  Structure and structure-preserving algorithms for plasma physics.  Phys. Plasmas 24, 055502\.
  * Morrison and Greene [1980] Morrison, P.J., Greene, J.M., 1980\.  Noncanonical Hamiltonian density formulation of hydrodynamics and ideal magnetohydrodynamics.  Phys. Rev. Lett. 45, 790–793.
  * Morrison and Hazeltine [1984] Morrison, P.J., Hazeltine, R.D., 1984\.  Hamiltonian formulation of reduced magnetohydrodynamics.  Phys. Fluids 27, 886–897.
  * Morrison and Updike [2023] Morrison, P.J., Updike, M.H., 2023\.  An inclusive curvature-like framework for describing dissipation: metriplectic 4-bracket dynamics.  arXiv:2306.06787 .
  * Nomizu [1971] Nomizu, K., 1971.  On the spaces of generalized curvature tensor fields and second fundamental forms .
  * Sato and Morrison [2024] Sato, N., Morrison, P.J., 2024\.  A collision operator for describing dissipation in noncanonical phase space.  arXiv , 2401.15086arXiv:cond-mat.stat-mech/2401.15086.
  * Sudarshan and Mukunda [1974] Sudarshan, E.C.G., Mukunda, N., 1974\.  Classical Dynamics: A Modern Perspective.  John Wiley and Sons.
  * Taylor [1992] Taylor, J.E., 1992.  II–Mean curvature and weighted mean curvature.  Acta Metall. Matter 40, 1475–1485.
  * Taylor and Cahn [1998] Taylor, J.E., Cahn, J.W., 1998\.  Diffuse interfaces with sharp comers and facets: Phase field models with strongly anisotropic surfaces.  Physica D 112, 381–411.
  * Thiffeault and Morrison [2000] Thiffeault, J.L., Morrison, P.J., 2000\.  Classification of Casimir invariants of Lie-Poisson brackets.  Physica D 136, 205–244.
  * Tran et al. [2024] Tran, B., Burby, J., Southworth, B., 2024.  On radiation hydrodynamics.  under preparation .
  * Zaidni et al. [2023] Zaidni, A., Benjelloun, S., Morrison, P.J., 2023.  On ideal and real dispersed two-phase flows: mixtures of particles, droplets, or bubbles in a fluid.  under preparation .