Next Article in Journal
Numerical Simulation of Entropy Generation with Thermal Radiation on MHD Carreau Nanofluid towards a Shrinking Sheet
Next Article in Special Issue
Chemical Reactions Using a Non-Equilibrium Wigner Function Approach
Previous Article in Journal
Distant Supervision for Relation Extraction with Ranking-Based Methods
Previous Article in Special Issue
Analysis of the Chaotic Behavior of the Lower Hybrid Wave Propagation in Magnetised Plasma by Hamiltonian Theory
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Hydrodynamic Theories for Flows of Active Liquid Crystals and the Generalized Onsager Principle

1
Beijing Computational Science Research Center, Beijing 100193, China
2
School of Mathematical Sciences and LPMC, Nankai University, Tianjin 300071, China
3
Departments of Mathematics and Biomedical Engineering, University of North Carolina at Chapel Hill, Chapel Hill, NC 27599, USA
4
Department of Mathematics, Interdisciplinary Mathematics Institute and NanoCenter at USC, University of South Carolina, Columbia, SC 29028, USA
*
Author to whom correspondence should be addressed.
These authors contributed equally to this work.
Entropy 2016, 18(6), 202; https://doi.org/10.3390/e18060202
Submission received: 5 December 2015 / Revised: 30 April 2016 / Accepted: 17 May 2016 / Published: 24 May 2016

Abstract

:
We articulate and apply the generalized Onsager principle to derive transport equations for active liquid crystals in a fixed domain as well as in a free surface domain adjacent to a passive fluid matrix. The Onsager principle ensures fundamental variational structure of the models as well as dissipative properties of the passive component in the models, irrespective of the choice of scale (kinetic to continuum) and of the physical potentials. Many popular models for passive and active liquid crystals in a fixed domain subject to consistent boundary conditions at solid walls, as well as active liquid crystals in a free surface domain with consistent transport equations along the free boundaries, can be systematically derived from the generalized Onsager principle. The dynamical boundary conditions are shown to reduce to the static boundary conditions for passive liquid crystals used previously.

1. Introduction

Active matter systems are abundant in Nature and man-made materials. They include many familiar systems like a bacterial suspension, a flock of birds, a school of fish, an ensemble of catalytic nanomotors, photo-powered nanoparticles, and the cytoskeleton in a live cell [1]. They share a common feature: the fundamental active “particles” are motile and can move on their own with an energy input created either internally (metabolic processes or ATP-motor activity) or externally from magnetic or electric fields, photonic, chemical gradients, or catalysis. In these diverse examples, active matter systems exhibit anomalous density fluctuations at scales far removed from single particles and emergent collective behavior. Dynamic phase transitions, spatio-temporal structures and symmetry-breaking induced dynamic patterns are prominent phenomena commonly seen in such systems [1,2,3]. Here, we are interested in flows of active liquid crystals, where the individual particles are self-propelling. We consider both polar and apolar self-propelling particles, where the passive equilibrium is nematic above a critical particle concentration. We refer to [1,2,3,4,5,6,7,8,9,10,11,12,13,14] for further details, and close the introduction with selected contributions of particular relevance to this paper. The review by Marchetti et al. [1] is especially recommended, from which we recall key elements and distinctions of active nematic “fluids”.
Models for “dry” active matter systems that do not conserve momentum have been explored by many, notably, by Chaté et al. [15,16,17,18,19] and Bertin et al. [20,21,22]. For wet active matter systems, a useful theoretical framework is a continuum model that generalizes liquid crystal hydrodynamics to include new features related to particle-scale activation [4]. For polar particles, continuum models employ a single vector to describe the broken head-tail symmetry; these polar vector models are identified as a low moment (zero and first) approximation of kinetic theory [1,2,23]. For apolar particles, models have employed either a single vector with imposed fore-aft symmetry (the nematic director) or a second order tensor (the nematic order tensor). For more general active liquid crystals whose active particles/molecules are primarily polar while the interactions may be dominated by apolar interactions, the first three moments (density, polarity vector, and the nematic order tensor) of the orientational distribution function are modeled either directly or by moment closures of the kinetic theory [3]. More generally, active liquid crystals are described by kinetic theories, generalizing the Doi-Edwards-Hess formalism for passive liquid crystals, cf. studies by Marchetti and Liverpool, Shelley, Saintillan et al. and Forest et al. [1,24,25,26,27,28,29,30]. In this formulation, microscopic symmetry is tacitly incorporated through the interaction potential, self-propelling velocity, as well as phenomenologically modeled active forces [31]. This formulation can unify all three types of continuum models, where any low moment model is recoverable from closure approximations. For further details, we refer to [1,2,4,32,33,34].
Joanny et al. derived a class of models for active polar liquid crystal gels using a formulation pioneered by Onsager for non-equilibrium thermodynamical systems [35]. In this paper, we extend the approach to more general hydrodynamic theories for active liquid crystal solutions and gels applicable to polar or apolar active matter systems [36]. We call this approach the generalized Onsager principle. First, we recall what is the Onsager theory for a matter system not far from equilibrium [37,38,39] and then define what we mean by the generalized Onsager principle.
We consider a matter system with a stable equilibrium at 0 , in which the fluctuations of a set of coarse-grained variables x = ( x 1 , , x n ) T are measured relative to their most probable (equilibrium) values x = 0 . The entropy of the system S, a function of x , reaches its maximum value S 0 at equilibrium. Expanding the entropy function in a Taylor series about the equilibrium, it can be approximated by the quadratic form
S = S 0 + Δ S ( x ) , Δ S = 1 2 x T · H · x ,
where H is the Hessian, a symmetric, positive-definite matrix when S has continuous second order derivatives at a maximum. The probability density at state x near the equilibrium is related to Δ S ( x ) via the Boltzmann distribution f ( x ) = f 0 e x p [ Δ S ( x ) / k B T ] , with k B the Boltzmann constant, T the absolute temperature and f 0 a normalization constant.
When the system deviates from equilibrium, spontaneous irreversible processes arise in response to the generalized thermodynamic force X conjugate to x :
X = ( Δ S x ) = H · x ,
which is linear in x due to the quadratic form of Δ S ( x ) in Equation (1). Thus, entropy production (or change of entropy) is given by
Δ S = 1 2 x T · X .
i.e., entropy production is determined by the inner product of the conjugate variable X and the original thermodynamic variable x .
For a small deviation from equilibrium, the system is assumed in the linear response regime so that the state x ( t ) evolves according to the kinetic equation
x ˙ = L · X ,
or, equivalently,
X = R · x ˙ ,
where the kinetic coefficients (mobility) L = ( L i j ) form a symmetric matrix and so do the friction coefficients R = ( R i j ) = L 1 , according to the Onsager reciprocal relation [37,40]. Off-diagonal entries of ( L i j ) and ( R i j ) are referred to as cross-coupling coefficients between different irreversible processes labeled by i and j. Under the condition that S is an even function of x , Onsager derived the reciprocal relation [37,38]
L i j = L j i ,
from the microscopic reversibility, that is for any t > 0 and τ,
x ( t ) x T ( t + τ ) = x ( t + τ ) x T ( t ) ,
where the ensemble average is taken with respect to the probability density function f ( x , t ) . We note that, for anisotropic and active matter systems, condition (7) may not be valid.
From the entropy function, we calculate the rate of entropy production as follows
S ˙ = x ˙ T · X = X T · L · X ,
which is given as an inner product of a generalized flux x ˙ and the generalized force X . In order for entropy to increase when approaching the stable equilibrium in an irreversible process, L must be nonnegative definite. We next summarize the three key ingredients in the Onsager theory for nonequilibrium systems near a stable equilibrium:
(1)
the kinetic Equation (4), which gives the governing system of equations for dynamics of the non-equilibrium system;
(2)
the Onsager reciprocal relation (6);
(3)
the nonnegative definiteness of L ensures entropy production while approaching the equilibrium in an irreversible process.
We term (1)–(3) the Onsager principle in place of the Onsager theory in this paper. This statement of the Onsager principle is equivalent to the Onsager maximum action principle for irreversible processes [37,38,39,40]. This is also equivalent to the minimum entropy production principle summarized by Prigogine for irreversible processes in stationary (or steady) states [41].
We next demonstrate that this is also equivalent to the statement on energy dissipation for nonequilibrium systems. Recall the first law of thermodynamics,
d E k + d U = d Q d W ,
where d E k is the increase of the kinetic energy of the matter system, d U is the increase of the internal energy of the matter system, d Q is the heat added to the system and d W is the work done to the surrounding by the system. We denote entropy S for the matter system as [42]
S = S i + S e ,
where S i is the entropy generated internally in the system and S e is the entropy exchanged with the surrounding. The Helmholtz free energy F of the system is defined by
F = U T S i ,
where T is the absolute temperature. So, it follows from Equation (9) that
d F = d U d T S i T d S i = d E k + d Q S i d T T d S i d W .
For an irreversible process,
T d S i = d Q + d Q ,
where d Q is the energy lost internally. Then,
d F + d E k = S i d T d W d Q .
For an isothermal system, d T = 0 . So,
d F + d E k = ( d W + d Q ) ,
This is known as the energy dissipation. Thus, the energy dissipation rate is given by
d d t ( E k + F ) = d ( U + E k ) d t T d S i d t = d W d t d Q d t = T X T · L · X .
For a closed system, d S e d t = 0 and d ( U + E k ) d t = 0 . So
d d t ( E k + F ) = T d S d t = d W d t d Q d t = T X T · L · X .
For an open system, d ( U + E k ) = T d S e , where it is assumed that there is no work done during the exchange. So,
d d t ( E k + F ) = T d S d t = d W d t d Q d t = T X T · L · X .
This establishes the relation between the energy dissipation rate and the entropy production rate, which applies to both the closed and open system.
In this paper, we generalize the Onsager principle to a non-equilibrium system in which the source of external or chemical energy that gives rise to self-propulsion in the active matter is accounted for. Specifically, we consider the dissipation rate defined by the combination of the kinetic energy and the Helmholtz free energy. In [35], this combination is called the free energy, which indeed serves as the role of free energy in the nonequilibrium systems. We remove the symmetry constraint imposed on the mobility coefficient in the Onsager principle and abandon the non-negativeness of the mobility coefficient imposed for irreversible processes to allow reversible systems. This extension is important for anisotropic viscoelastic material systems including active matter systems, where both irreversible and reversible process can coexist. For these material systems, the mobility coefficient L can be decomposed into a symmetric and an antisymmetric part. The antisymmetric component, representing reversible dynamics, does not contribute to energy dissipation while the symmetric part does. The symmetric part, representing irreversible dynamics, is required to be nonnegative definite to ensure energy dissipation. We remark that the antisymmetric part corresponds to the Casimir principle while the symmetric part corresponds to the Onsager principle [42,43].
In the active matter system, we consider the source of chemical energy (e.g., the one produced by ATP hydrolysis in myosin motors or through catalytic reaction) or external source converted chemical energy (e.g., energy due to photo induced chemical reactions or magnetic fields) as a source for an additional free energy. We first derive the governing system of equations in a domain with a fixed boundary together with the consistent boundary conditions that do not contribute to energy dissipation. Then, we extend it to a two-phase case in which the active liquid crystal system resides in a domain separated from a passive fluid matrix by free boundaries. In this case, we also consider the surface transport phenomena which may exist in such a matter system. In both cases, we let the active component of the free energy enters into the total stress through a reversible process and the dynamics of the internal variables (density, polarity vector, and nematic tensor) through irreversible processes. Then, the resultant theories of active matter systems may not respect energy dissipation anymore but their corresponding passive component, the limit when the activity is absent, must. This is the basic principle that we follow in the following derivation.
These models of complex active matter systems have applications to cell motility, a school of fish, a colony of bacteria in a viscous fluid bath, and the cortical layer in a live cell, guaranteeing that the system conserves momentum and the corresponding passive component satisfies the correct energy dissipation principle. In this way, the detailed physical and chemical processes of self-propulsion are determined specific to the active matter system and incorporated in the framework of the generalized Onsager principle. We proceed in this paper to illustrate this framework for active liquid crystals. We remark that there exist other approaches for developing extended hydrodynamic theories for complex fluid systems like the generalized Poisson bracket approach discussed in the book [44], the GENERIC formalism advocated in [45,46], and several more nonequilibrium thermodynamic frameworks reviewed by Sagis and others in [47,48,49,50]. All these approaches agree in principle on the fundamental mathematical structure of the non-equilibrium thermodynamic theory. In particular, the approach that has been used by Sagis et al. is especially close to what we are using in this paper [47,48,49]. However, our approach presented in this paper is simpler, constructive, and straight-forward in that it yields (derives or stipulates) the constitutive relations explicitly that warrants energy dissipation for the passive components for both the reversible and the irreversible nonequilibrium processes.
The rest of the paper is organized as follows. In Section 2, we derive the hydrodynamic model along with the boundary conditions in a fixed domain for an active liquid crystal system, closing the model on the first three moments of a probability distribution function. In Section 3, we present the model for an apolar active liquid crystal system. In Section 4, we extend the study to domains with free surfaces. Then we provide a closing remark.

2. A General Hydrodynamic Model for Active Liquid Crystals

We systemically derive a general hydrodynamic model for flows of active liquid crystals using the generalized Onsager principle at the continuum scale. We assume that the active liquid crystal system is composed of two components: the first component consists of the active anisotropic particles with a mass density ρ 1 and velocity v 1 and the second component consists of the solvent which has a mass density ρ 2 with velocity v 2 . The corresponding mass conservation laws for the two components are given respectively by
t ρ 1 + · ( ρ 1 v 1 ) = 0 , t ρ 2 + · ( ρ 2 v 2 ) = 0 .
Introducing the total mass density ρ = ρ 1 + ρ 2 and the mass-averaged center-of-mass velocity v = ( ρ 1 v 1 + ρ 2 v 2 ) / ρ , the total linear momentum density (or density flux) can be expressed as ρ v . The total mass is conserved:
t ρ + · ( ρ v ) = 0 .
For the active particle component, we rewrite the current into a convective part moving with the center-of-mass velocity v and a diffusive part associated to the relative flux between the two components:
t ρ 1 + · ( ρ 1 v + j ^ ) = 0 ,
where j ^ = ρ 1 ρ 2 ( v 1 v 2 ) ρ is the diffusive current. In the following, we define the active particle number density as c = ρ 1 m 1 , where m 1 is the mass of the active particle. We denote its flux as j = j ^ m 1 . The transport equation for c can then be written as
t c + · ( c v + j ) = 0 .
To describe the system’s internal properties (or microstructure), a hierarchy of order parameters is necessary. A scalar concentration variable is needed to describe the density fluctuation in space and time. A vector order parameter is needed to describe the average polarity. A rank 2 tensor order parameter is needed to characterize the nematic order or the orientational correlation among the particles. So, the minimal number of order parameters in any viable model should be three: a scalar, a vector, and a second order tensor. We denote the polar order by a vector field p and the nematic order by a second order tensor field Q , in addition to the density variable c.
Consider the active particle system as a discrete system with N particles located at r n ( t ) , n = 1 , , N , respectively, and each of the particles has a self-propelled velocity ν ^ n ( t ) , n = 1 , , N . We can normalize ν ^ n ( t ) as a unit vector to represent the direction of motion and describe speed fluctuations with a scalar, but we will not employ this normalization in this paper. We assume the domain of ν ^ , V , is a compact domain in R 3 . The one-particle phase-space distribution f d ( r , ν ^ , t ) is defined as follows:
f d ( r , ν ^ , t ) = n δ ( r r n ( t ) ) δ ( ν ^ ν ^ n ( t ) ) .
The number density of active particles c ( r , t ) is the zeroth moment
c ( r , t ) = V f d ( r , ν ^ , t ) d ν ^ = n δ ( r r n ( t ) ) .
The first moment of the distribution function is given by
c ( r , t ) p ( r , t ) = V ν ^ f d ( r , ν ^ , t ) d ν ^ = n ν ^ n ( t ) δ ( r r n ( t ) ) .
The polarization vector field p ( r , t ) follows as a normalized first moment. The second moment is defined by
c ( r , t ) Q ( r , t ) = V ( ν ^ ν ^ ν 2 d ) f d ( r , ν ^ , t ) d ν ^ = n ( ν ^ n ( t ) ν ^ n ( t ) ν ^ n 2 d ) δ ( r r n ( t ) ) ,
where d is the system’s dimensionality, from which a nematic order tensor Q ( r , t ) can be defined as a normalized second moment. Notice that the tensor Q is traceless.
Often, one prefers the moments defined as primitive variables for developing multiphase models. In this context, the polarity vector and the nematic order tensor can then be defined as follows:
p ( r , t ) = V ν ^ f d ( r , ν ^ , t ) d ν ^ = n ν ^ n ( t ) δ ( r r n ( t ) ) , Q ( r , t ) = V ( ν ^ ν ^ ν 2 d ) f d ( r , ν ^ , t ) d ν ^ = n ( ν ^ n ( t ) ν ^ n ( t ) ν ^ n 2 d ) δ ( r r n ( t ) ) .
The free energy of the active particle system in volume element V must be prescribed in terms of the internal variables, in this case the first three moments of f d and their gradients (i.e., we assume weak non-locality for this paper):
F = F ( c , c , p , p , Q , Q ) = V f ( c , c , p , p , Q , Q ) d x ,
where f is the free energy density. At the boundary of the fixed material domain V , we have an interfacial free energy
G = G ( c , p , Q ) = V g ( c , p , Q ) d S .
The active free energy in the active matter system is denoted as
A = V R ( x , t ) μ ^ d x ,
where R ( x , t ) is the number density of the fundamental energy generating units in the active matter system and μ ^ is the energy gain per fundamental unit at position x , time t. For instance, μ ^ is the energy gain per ATP molecule for F-actin filaments or microtubule in a live cell while R ( x , t ) is the number density of the ATP molecules. The sum of the kinetic and free energy of the system is called the total energy here or extended free energy in [35] given by
E t o t a l = V [ ρ 2 v 2 ] d x + F + G + A .
Since we are interested in an isothermal system, the temperature is assumed a constant and so the energy dissipation rate at a constant absolute temperature T is given by
d E t o t a l d t = V d x { 1 2 ( ρ v 2 ) t + f t r μ ^ } + V g t d S .
There are four parts in the integration: the first and second parts are the rate changes of the kinetic energy and the active particle free energy, respectively, the third one is the energy reduction rate of the chemical energy, where r = t R is the consumption rate (the number of the fundamental energy producing units consumed per unit time and unit volume), the fourth is the rate of change of the surface free energy. The differential of the free energy density is given by
f t = f c c t + f p · p t + f Q : Q t + f ( c ) · ( c ) t + f ( p ) : ( p ) t + f ( Q ) ( Q ) t = μ c t h · p t G : Q t + · ( f ( c ) c t + f ( p ) · p t + f ( Q ) : Q t ) .
Now we identify the conjugate variables to c , p , and Q through variations of the free energy density f. The field conjugate to the density c is the chemical potential μ = δ f δ c = f c · f ( c ) , the field conjugate to the polarization p is the molecular field h = δ f δ p = f p + · f ( p ) , and the field conjugate to the nematic order Q is the variation G = δ f δ Q = f Q + · f ( Q ) , another molecular field. If Q is traceless, G = [ 1 2 ( δ f δ Q + δ f δ Q T ) I 3 t r ( δ f δ Q ) ] . We denote the product of two second order tensors as A : B = A α β B α β . So
d E t o t a l d t = V d x { 1 2 ( ρ v 2 ) t + μ c t h · p t G : Q t r μ ^ } + V d S { ( f ( c ) c t + f ( p ) · p t + f ( Q ) : Q t ) · n } + V d S { ( g c c t + g p · p t + g Q : Q t ) } ,
where V is the surface of volume V and n is the external unit normal of V . We assume the following so that the surface integral is zero, i.e., the surface does not contribute to energy dissipation of the system:
n · f c + g c = 0 , n · f p + g p = 0 , n · f Q + g Q = 0 .
These define the boundary conditions of the active matter system at the solid walls for the three internal variables c , p and Q . From the previous discussion, the conservation laws of the active particle number, total system mass and system momentum are given by
c t + · ( c v + j ) = 0 , ρ t + · ( ρ v ) = 0 , ( ρ v ) t = · ( σ ρ vv ) ,
where j is the diffusive current of the active particles, σ is the momentum flux. Then, we have
V d x { 1 2 ( ρ v 2 ) t + μ c t } = V d x { β v α σ α β + c v α α μ + j α α μ } + V d S { v · ( 1 2 ρ v 2 I c μ I + σ ) · n μ j · n } .
If we assume v = 0 and j · n = 0 at the boundary, the last surface integration is zero. Once again, the surface does not contribute to energy dissipation of the active matter system. These define additional boundary conditions for the velocity and the internal variables via the excessive flux. If the surface term is not assigned into zero, it will contribute to the energy dissipation. The latter case includes moving boundaries, which we will not pursue in this study.
We denote
D α β = 1 2 ( α v β + β v α ) , Ω α β = 1 2 ( α v β β v α ) , p ˙ = t p + v · p + Ω · p , Q = t Q + v · Q + [ Ω · Q Q · Ω ] ,
where D is the strain rate tensor, Ω is the vorticity tensor of the velocity field v , p ˙ is the convected co-rotational derivative of vector p , Q is the co-rotational derivative of tensor Q . Then
h · p t = h α p ˙ α + v α h γ α p γ + 1 2 β v α ( p α h β p β h α ) , G : Q t = G α β Q α β + v α G β γ α Q β γ + β v α ( Q α γ G γ β G α γ Q γ β ) ,
where G : [ Ω · Q Q · Ω ] = β v α ( Q α γ G γ β G α γ Q γ β ) due to symmetry of Q , G and antisymmetry of Ω. We introduce the antisymmetric part of the stress σ a and the Ericksen stress σ e (Refer to Appendix A) respectively as follows:
σ α β a = 1 2 ( p α h β p β h α ) + ( Q α γ G γ β G α γ Q γ β ) , σ α β e = ( f c μ ) δ α β f ( β c ) α c f ( β p γ ) α p γ f ( β Q γ θ ) α Q γ θ , β σ α β e = ( c α μ + h γ α p γ + G β γ α Q β γ ) .
Then the energy dissipation is simplified into
d E t o t a l d t = V d x { β v α [ σ α β σ α β e σ α β a ] + j α α μ h α p ˙ α G α β Q α β r μ ^ } V d x β ( v α σ α β e ) .
From the vanishing surface terms, we summarize the boundary conditions as follows
v | V = 0 , f ( c ) · n + g c = 0 , j · n = 0 , f ( p ) · n + g p = 0 , f ( Q ) · n + g Q = 0 .
For c , p , Q , the boundary conditions (BCs) are mixed (Robin). There are two boundary conditions for c due to the nature of the transport (the Cahn-Hilliard equation) of c. We remark that if we do not postulate the transport equation of c in a conservative form (e.g., we use an Allen–Cahn type equation to describe the transport of c), the boundary condition on j will not be necessary. If the surface energy density dominates, we arrive at the Dirichlet BCs (or the anchoring boundary condition) for the internal variables:
g c = 0 , g p = 0 , g Q = 0 .
There are two issues to be addressed in order to complete the derivation through the generalized Onsager principle. One is to derive the Gibbs–Duhem relation to identify the Ericksen stress σ α β e . The other one is to prove that the stress σ s , given by
σ α β s = σ α β σ α β e σ α β a ,
is in fact symmetric. We address these two issues in Appendix A to keep the presentation flowing. For the time being, we assume both are true. Then, we rewrite the energy dissipation functional as follows
d E t o t a l d t = V d x { D α β σ α β s + j α ( α μ ) + h α p ˙ α + G α β Q α β + r μ ^ } .
The generalized fluxes and the corresponding generalized forces, identified from the energy dissipation, are listed below
Flux Force σ α β s D α β Q α β G α β p ˙ α h α j α α μ r μ ^ .
Note that, the generalized forces have different signatures under time inversion. While D changes sign, the other forces do not. Therefore we distinguish between the components of the fluxes that show the same behavior under time inversion as the dissipative conjugated forces or fluxes, and those that show the opposite behavior called reactive conjugated forces or fluxes. We denote the various components by superscripts d and r, respectively.
The phenomenological equations for the dissipative currents can be written as
σ α β s , d Q α β d p ˙ α d j α d r d = A 0 0 0 0 0 0 1 γ 2 δ α k δ β l χ 1 2 A α β k χ 2 2 A α β k A 5 , α β 0 χ 1 2 A α k l T 1 γ 1 δ α k λ δ α k B α 0 χ 2 2 A α k l T λ δ α k γ δ α k ( κ + ω 0 c ) p α 0 A 5 , k l B k ( κ + ω 0 c ) p k Λ D k l G k l h k k μ μ ^ ,
where
A α β k = ( p α δ β k + p β δ α k ) 2 3 δ α β p l δ l k , A α k l T = 2 ( p β δ α k δ β l p β 3 δ α β δ k l ) , A 0 = 2 η δ α k δ β l + ( η ¯ 2 3 η ) δ α β δ k l + α 1 ( Q α k δ β l + δ α k Q β l ) + α 2 Q k l Q α β + α 3 ( p α p k δ β l + δ α k p β p l ) + α 4 p k p l p α p β , B α = λ 1 p α ω p · p α ω p α · p , A 5 = ζ 1 ( pp p 2 3 I ) + ( ζ 2 ω 2 · p ω 2 p · ) Q .
It is easy to show that the coefficient matrix is symmetric. We insist the passive limit of the active matter system is a dissipative system. Then, the left upper submatrix is nonnegative definite. So, the diagonal parameters η , γ , γ 1 , γ 2 are nonnegative together with a set of constraints on the model parameters independent of μ ^ . The reactive terms can be written as
σ α β s , r Q α β r p ˙ α r j α r r r = 0 A 1 A 2 A 3 A 4 A 1 0 0 0 0 A 2 0 0 0 0 A 3 0 0 0 0 A 4 0 0 0 0 D k l G k l h k k μ μ ^ ,
where
A 1 = ν 0 δ α k δ β l + a [ Q α k δ β l + δ α k Q β l ] a ( Q k l ( Q α β + 1 3 δ α β ) ) + θ 1 δ k l δ α β , A 2 = ν 1 2 ( p β δ α k + p α δ β k ) + θ 2 p k δ α β , A 2 = ν 1 p β δ α k δ β l + θ 2 p α δ k l , A 3 = ν 2 2 ( p β δ α k + p α δ β k ) + θ 3 p k δ α β , A 3 = ν 2 p β δ α k δ β l + θ 3 p α δ k l , A 4 = ζ p α p β + ζ 0 Q α β + β 1 ( α p β + β p α ) + β 2 p · Q α β + β 3 2 Q α β ( ζ 2 + ζ 3 p α p α ) δ α β , A 4 = A 4 δ α k δ β l .
For the reactive fluxes, the coefficient matrix is antisymmetric, implying Flux r · Force = 0 ; so, the reactive processes do not contribute to energy dissipation.
We summarize the equations in the following
t c + · ( ( v + ω 0 μ ^ p ) c ) = · ( γ μ + χ 2 G · p χ 2 3 t r ( G ) p + λ h + κ p μ ^ + ν 2 D · p + θ 3 D γ γ p ) , t p + · ( ( v + ω μ ^ p ) p ) · v p + Ω · p = 1 γ 1 h + χ 1 G · p χ 1 3 t r ( G ) p λ μ + λ 1 μ ^ p + ν 1 D · p + θ 2 D γ γ p , Q t + · ( ( v + ω 2 μ ^ p ) Q ) · v Q + Ω · Q Q · Ω a [ Q · D + D · Q ] = 1 γ 2 G + χ 1 2 ( ph + hp ) χ 2 2 ( p μ + μ p ) + ζ 1 μ ^ pp + ν 0 D a ( Q : D ) Q + ( θ 1 D γ γ a Q : D 3 χ 1 p · h 3 + χ 2 p · μ 3 μ ^ ζ 1 p 2 3 ) I + ζ 2 μ ^ Q , t ρ + · ( ρ v ) = 0 , t ( ρ v ) + · ( ρ vv ) = · σ , σ = σ s + σ e + σ a , σ α β e = ( f c μ ) δ α β f ( β c ) α c f ( β p γ ) α p γ f ( β Q γ θ ) α Q γ θ , · σ e = ( c μ + p · h + Q : G ) , σ a = 1 2 ( ph hp ) + ( Q · G G · Q ) , σ s = Π I + 2 η D + α 1 ( Q · D + D · Q ) + α 2 ( Q : D ) Q + α 3 ( pp · D + D · pp ) + α 4 ( pp : D ) pp a ( Q · G + G · Q ) ν 0 G + a ( Q : G ) Q ν 1 2 ( ph + hp ) + ν 2 2 ( p μ + μ p ) + ζ μ ^ pp + ζ 0 μ ^ Q + β 1 μ ^ ( p + p T ) + β 2 μ ^ p · Q + β 3 μ ^ 2 Q ,
where Π = ( η ¯ 2 3 η ) D γ γ a Q k l G k l 3 + θ 1 G γ γ + θ 2 p γ h γ θ 3 p γ γ μ + ( ζ 2 + ζ 3 p γ p γ ) μ ^ is the pressure. If we assume the total density ρ = ρ 0 is constant, that means the fluid is incompressible, · v = 0 , D γ γ = 0 . Then we can drop the parameters θ i , i = 1 , 2 , 3 and Π is the hydrostatic pressure. The boundary conditions are given by Equation (38).
All the terms related to μ ^ are active energy contributions in the system. If the active energy is removed from the system, the system is reduced to the traditional passive liquid crystal system, and our model recovers the usual liquid crystal hydrodynamics model [5,51,52].

3. Hydrodynamic Model for Apolar Active Liquid Crystal Systems

For the apolar system, we assume the polar vector p = 0 and only consider the tensor order parameter Q and concentration c. Then, the model reduces to the following.
t c + · ( v c ) = · ( γ μ ) , t Q + v · Q + Ω · Q Q · Ω a [ Q · D + D · Q ] = 1 γ 2 G + ν 0 D a ( Q : D ) ( Q + I 3 ) + ζ 2 μ ^ Q , ρ t + · ( ρ v ) = 0 , t ( ρ v ) + · ( ρ vv ) = · ( σ s + σ e + σ a ) , σ α β e = ( f c μ ) δ α β f ( β c ) α c f ( β Q γ θ ) α Q γ θ , · σ e = ( c μ + Q : G ) , σ a = ( Q · G G · Q ) , σ s = Π I + 2 η D + α 1 ( Q · D + D · Q ) + α 2 ( Q : D ) Q a ( Q · G + G · Q ) ν 0 G + a ( Q : G ) Q + ζ 0 μ ^ Q + β 3 μ ^ 2 Q .
where Π = a Q k l G k l 3 + ζ 2 μ ^ .
The free energy functional is given by [1]
F = r { α ( c ) 2 Q : Q + β ( c ) 3 t r ( Q 3 ) + γ ( c ) 4 ( Q : Q ) 2 + K ( c ) 2 ( Q Q ) + L ( c ) 2 · Q 2 + C ( c ) Q : δ c c 0 + A ( c ) 2 ( δ c c 0 ) 2 + B ( c ) 2 c 2 } .
The system is nematic at equilibrium if α < 0 and γ > 0 . K is the Frank elasticity coefficient, the C term is the bilinear coupling of Q and c. A is the compression modulus of density fluctuations δ c = c c 0 . We can generalize A ( c ) 2 δ c 2 into a function of c: h ( c ) .
The conjugate fields are given by
μ = δ F δ c = α ( c ) 2 Q : Q + β ( c ) 3 t r ( Q 3 ) + γ ( c ) 4 ( Q : Q ) 2 + K ( c ) 2 ( Q Q ) + L ( c ) 2 · Q 2 + C ( c ) Q : δ c c 0 + A ( c ) 2 ( δ c c 0 ) 2 + B ( c ) 2 c 2 + C c 0 : Q + A δ c c 0 2 · ( B c ) , G α β = δ F δ Q α β = [ α Q + β ( Q 2 Q : Q I 3 ) + γ Q : Q Q · ( K Q ) ( 1 2 α ( L i Q i β ) + 1 2 β ( L i Q i α ) α ( L i Q i α ) 3 I ) + C c 0 ( c t r ( c ) I 3 ) ] .

4. Hydrodynamic Model for Free Surface Flows of Active Matter Systems

4.1. Models for Passive Liquid Crystals

We first consider a mixture of passive liquid crystals in which a free surface separates the liquid crystal from the other fluids. We denote the total energy of the liquid crystal system in a given material domain V by
E = V [ ρ 2 v 2 + f ( c , c , p , p , Q , Q ) ] d x + Γ [ ρ s 2 v 2 + g ( c s , s c s , p s , s p s , Q s , s Q s ) ] d s ,
where f is the bulk free energy density, g is the “excessive” [53,54] surface free energy density at a material interface Γ within domain V, defined by Φ ( x , t ) = 0 , ρ s is the surface excessive mass density, c s is the surface excessive active liquid crystal concentration, p s and Q s are the polarity vector and nematic tensor on the surface, and s = I s · is the surface gradient operator, where I s = I nn , n is the unit normal of the interfacial surface, given by n = Φ Φ . Here, we assume the excessive free energy density g on the free surface depends on the gradients of the order parameters as well. We assume the velocity across the interface is continuous so that
( t + v · ) Φ = 0 .
The bulk conservation laws are the mass, active species and momentum conservation:
t ρ + · ( ρ v ) = 0 , t ( ρ v ) = · ( σ ρ vv ) + F e , t c + · ( c v + j ) = 0 ,
where σ is the bulk stress tensor, F e is the bulk elastic force and j is the diffusive flux for the liquid crystal.
The surface conservation laws are surface excessive mass, active species conservation and surface momentum conservation:
t ρ s + s · ( ρ s v ) = 0 , t ( ρ s v ) = s · ( σ s ρ s vv ) + F s , t c s + s · ( c s v + j s ) = j c ,
where σ s is the excessive surface stress, F s is the excessive surface elastic force and j s is the excessive surface diffusive flux for liquid crystals. We note that v ( x , t ) | Γ = v ( x s , t ) because the surface is a material surface. We assume v | V = 0 in the following derivation and define
Ω s = 1 2 [ s v · I s I s · s v T ] , D s = 1 2 [ s v · I s + I s · ( s v ) T ] .
We define the invariant derivatives on the surface as follows
d s d t c s = t c s + v · s c s , p ˜ s = t p s + v · s p s + Ω s · p s , Q ^ s = t Q s + v · s Q s + Ω s · Q s · I s I s · Q s · Ω s .
In the following, we give a simple proof for the invariant derivative: p ˜ s . If both the background moving surface and the director p s on the surface rotate around the normal n at the same angular velocity θ. Then the surface velocity v s = I s · ( θ × r ) , and Ω s · p s = I s · ( θ × p s ) . At time t, the director at location r is p s ( r , t ) . After a small time Δ t , at time t + Δ t ,the director at the location r + I s · ( θ × r ) Δ t is p s ( r + I s · ( θ × r ) Δ t , t + Δ t ) = p s ( r , t ) + I s · ( θ × p s ( r , t ) ) Δ t + O ( Δ t 2 ) ,
p s t = lim Δ t 0 p s ( r , t + Δ t ) p s ( r , t ) Δ t = lim Δ t 0 p s ( r , t + Δ t ) p s ( r + I s · ( θ × r ) Δ t , t + Δ t ) + p s ( r + I s · ( θ × r ) Δ t , t + Δ t ) p s ( r , t ) Δ t = lim Δ t 0 I s · ( θ × r ) · p s ( r + I s · ( θ × r ) a Δ t , t + Δ t ) Δ t + I s · ( θ × p s ( r , t ) ) Δ t Δ t , a ( 0 , 1 ) = I s · ( θ × r ) · p s ( r , t ) + I s · ( θ × p s ) = v · s p s Ω s · p s .
So, p ˜ s = t p s + v · s p s + Ω s · p s = 0 under a pure rotation on the tangent plane, which is equivalent to say that p ˜ s is an invariant derivative. For tensor Q s , we can show its rotational invariance in the tangent plane analogously using its definition.
Then the energy dissipation rate is calculated as follows.
d E d t = V [ ( v ) T : ( σ ) + v · F e + δ f δ c c t + δ f δ p · p t + δ f δ Q : Q t + · ( v f ) ] d x + Γ [ [ f i p j t p j n i + f i Q j k t Q j k n i + f i c t c n i + vn : σ ] ] d s + Γ [ g c s t c s + g p s · t p s + g Q s : t Q s + g s c s · t ( s c s ) + g s p s : t ( s p s ) + g s Q s t ( s Q s ) + s · ( v g ) ] d s + Γ [ v · ( s · σ s + F s ) ] d s = V [ ( v ) T : ( σ σ a σ e ) + v · F e + j · μ h · p ˙ G : Q ] d x + Γ [ [ f i p j t p j n i + f i Q j k t Q j k n i + f i c t c n i + vn : σ μ ( c v i + j i ) n i v i σ i j e n j + f v i n i ] ] d s + Γ [ g c s d s c s d t + g p s · d s p s d t + g Q s : d s Q s d t + g s c s · t ( s c s ) + g s p s : t ( s p s ) + g s Q s t ( s Q s ) + g s c s · v · s ( s c s ) + g s p s : ( v · s ) ( s p s ) + g s Q s ( v · s ) ( s Q s ) + g s · v + s · ( v · σ s ) ] d s + Γ [ ( s v ) T : ( σ s ) + v · F s ] d s .
Here, the notation [ [ · ] ] denotes the jump of the variable { · } across the free surface. For simplicity, we assume the surface energy is short-ranged: g = g ( c s , p s , Q s ) ,
d E d t = V [ ( v ) T : ( σ σ a σ e ) + v · F e + j · μ h · p ˙ G : Q ] d x + Γ [ [ f i c t c n i + f i p j t p j n i + f i Q j k t Q j k n i + vn : σ μ ( c v i + j i ) n i v i σ i j e n j + f v i n i ] ] d s + Γ [ g c s d s c s d t + g p s · p ˜ s + g Q s : Q ^ s g p s · ( Ω s · p s ) g Q s : ( Ω s · Q s · I s I s · Q s · Ω s ) + g s · v + s · ( v · σ s ) ] d s + Γ [ ( s v ) T : ( σ s ) + v · F s ] d s .
We note that the one-sided limit of the internal variables, assuming they are defined globally, can be related to the surface value via the adsorption coefficients as follows:
c ± | Γ = K c , ± c s , p ± | Γ = K p , ± p s , Q ± | Γ = K Q , ± Q s ,
where K i , ± , i = c , p , Q are the adsorption coefficients from the + and − side of the surface. In the following, we consider one side ( ) of the surface as an isotropic viscous fluid where K i , = 0 , we only need one absorption coefficient, that is K i , + = K i . Hence, we drop the subscript ± on the K s . These define the relations between the internal variables on the surface and the ones in the bulk. With these assumptions, the energy dissipation reduces to
d E d t = V [ ( v ) T : ( σ σ a σ e ) + v · F e + j · μ h · p ˙ G : Q ] d x + Γ [ [ vn : ( σ σ e ) μ ( c v + j ) · n n · f p · ( v · s ) p n · f Q : ( v · s ) Q n · f c ( v · s ) c + f v · n ] ] d s + Γ [ ( g c s + K c n · f c ) d s c s d t + ( g p s + K p n · f p ) · p ˜ s + ( g Q s + K Q n · f Q ) : Q ^ s ( g p s + K p n · f p ) · ( Ω s · p s ) ( g Q s + K Q n · f Q ) : ( Ω s · Q s · I s I s · Q s · Ω s ) + g s · v + s · ( v · σ s ) ] d s + Γ [ ( s v ) T : ( σ s ) + v · F s ] d s .
It follows from the transport equation for c s that d s c s d t = s · v c s s · j s + j c , where j c is an excessive surface density growth rate. The energy dissipation reduces to
d E d t = V [ ( v ) T : ( σ σ a σ e ) + v · F e + j · μ h · p ˙ G : Q ] d x + Γ [ [ vn : ( σ σ e ) μ ( c v + j ) · n n · f p · ( v · s ) p n · f Q : ( v · s ) Q n · f c ( v · s ) c + f v · n ] ] d s + Γ [ ( g c s + K c n · f c ) ( s · v c s + s · ( I s · j s ) j c ) + ( g p s + K p n · f p ) · p ˜ s + ( g Q s + K Q n · f Q ) : Q ^ s ( g p s + K p n · f p ) · ( Ω s · p s ) ( g Q s + K Q n · f Q ) : ( Ω s · Q s · I s I s · Q s · Ω s ) + g s · v + s · ( v · σ s ) ] d s + Γ [ ( s v ) T : ( σ s ) + v · F s ] d s ,
where we have used the assumption j s = I s · j s . The integral term consists of
Γ s · ( v · σ s ) d s = Γ ( s · ( I s · v · σ s ) 2 κ v · σ s · n ) d s = Γ v · σ s · dL + Γ 2 κ v · σ s · n d s ,
where κ = 1 2 ( κ 1 + κ 2 ) is the mean surface curvature, κ 1 , κ 2 are the principal curvatures of the surface. The first term is assumed 0 because if Γ V it is true automatically, otherwise, Γ is a closed surface so that Γ = { ϕ } . Then, the second term is added to F s .
From Equation (62), we define μ s = g c s + K c n · f c , h s = ( g p s + K p n · f p ) , G s = ( g Q s + K Q n · f Q ) , which are the chemical potential, molecular field and the corresponding variation field to Q s on the surface. We then define the following
σ s = σ σ a σ e , σ s s = σ s ( g c s μ s ) I s 1 2 I s · ( p s h s h s p s ) ( I s · Q s · I s · G s I s · G s · I s · Q s ) , F s + [ [ ( σ σ e ) ] ] · n + [ [ f μ c ] ] n K c ( n · f c ) s c s K p s p s · ( n · f p ) K Q s Q s : ( n · f Q ) 2 κ σ s · n = 0 , [ [ μ j ] ] · n + j c μ s 0 , F e = 0 .
The energy dissipation reduces to
d E d t V [ ( v ) T : σ s + j · μ h · p ˙ G : Q ] d x + Γ [ ( s v ) T : σ s s + j s · s μ s h s · p ˜ s G s : Q ^ s ] d s .
Here we drop the term Γ s · ( I s · j s μ s ) d s , which is a line integral assumed to be 0. This formulation allows us to apply the generalized Onsager principle to obtain the constitutive relations. Since the passive liquid crystal system is a special limit of the more general active liquid crystal system, we will focus on the constitutive equations for the active liquid crystal.

4.2. Models for Active Liquid Crystals

For active liquid crystals, we need to add the active energy contribution r μ ^ to the total free energy functional. Namely,
d E d t V [ ( v ) T : σ s + j · μ h · p ˙ G : Q r μ ^ ] d x + Γ [ ( s v ) T : σ s s + j s · s μ s h s · p ˜ s G s : Q ^ s r s μ s ^ ] d s ,
where r is the active energy consumption rate and μ ^ is the active energy potential in the bulk, r s and μ s ^ are the corresponding quantities on the moving interface. In the bulk, we have proved that σ s is symmetry (see Appendix A), so that
( v ) T : σ s = D : σ s .
On the surface, we assume that σ s s = I s · σ s s · I s and it is symmetry, such that
( s v ) T : σ s s = D s : σ s s .
Then, we can apply the generalized Onsager principle to arrive at the constitutive relations. We denote the generalized flux as F = ( σ s , Q , p ˙ , j , r , σ s s , Q ^ s , p ˜ s , j s , r s ) and the generalized force as U = ( D , G , h , μ , μ ^ , D s , G s , h s , s μ s , μ s ^ ) . We find that the bulk free energy also gives a contribution to the surface generalized forces in G s , h s , s μ s . The generalized Onsager principle states that
F = ( D s y m + D a n t i ) · U ,
where
D s y m = D s y m b 0 0 D s y m s , D a n t i = D a n t i b 0 0 D a n t i s ,
D s y m b , D s y m s are 5 × 5 symmetric matrices and D a n t i b , D a n t i s are 5 × 5 antisymmetric matrices, respectively. We use subscripts α β for tensors and α for vectors in F and use subscripts k l for tensors and k for vectors in U . Then, the corresponding coefficients can be written as follows. First, the coefficients in the bulk are given by
D s y m b = A 0 0 0 0 0 0 1 γ 2 δ α k δ β l χ 1 2 A α β k χ 2 2 A α β k A 5 , α β 0 χ 1 2 A α k l T 1 γ 1 δ α k λ δ α k B α 0 χ 2 2 A α k l T λ δ α k γ δ α k ( κ + ω 0 c ) p α 0 A 5 , k l B k ( κ + ω 0 c ) p k Λ , D a n t i b = 0 A 1 A 2 A 3 A 4 A 1 0 0 0 0 A 2 0 0 0 0 A 3 0 0 0 0 A 4 0 0 0 0 ,
where
A α β k = ( p α δ β k + p β δ α k ) 2 3 δ α β p l δ l k , A α k l T = 2 ( p β δ α k δ β l p β 3 δ α β δ k l ) , A 0 = 2 η δ α k δ β l + ( η ¯ 2 3 η ) δ α β δ k l + α 1 ( Q α k δ β l + δ α k Q β l ) + α 2 Q k l Q α β + α 3 ( p α p k δ β l + δ α k p β p l ) + α 4 p k p l p α p β , B α = λ 1 p α ω p · p α ω p α · p , A 1 = ν 0 δ α k δ β l + a [ Q α k δ β l + δ α k Q β l ] a ( Q k l ( Q α β + 1 3 δ α β ) ) + θ 1 δ k l δ α β , A 2 = ν 1 2 ( p β δ α k + p α δ β k ) + θ 2 p k δ α β , A 2 = ν 1 p β δ α k δ β l + θ 2 p α δ k l , A 3 = ν 2 2 ( p β δ α k + p α δ β k ) + θ 3 p k δ α β , A 3 = ν 2 p β δ α k δ β l + θ 3 p α δ k l , A 4 = ζ p α p β + ζ 0 Q α β + β 1 ( α p β + β p α ) + β 2 p · Q α β + β 3 2 Q α β ( ζ 2 + ζ 3 p α p α ) δ α β , A 4 = A 4 δ α k δ β l , A 5 = ζ 1 ( pp p 2 3 I ) + ( ζ 2 ω 2 · p ω 2 p · ) Q .
Similarly, the coefficients of the matrices defined on the surface terms are given by
D s y m s = A 0 s 0 0 0 0 0 1 γ 2 s δ α k δ β l χ 1 s 2 A α β k s χ 2 s 2 A α β k s A 5 , α β s 0 χ 1 s 2 A α k l s , T 1 γ 1 s δ α k λ s δ α k B α s 0 χ 2 s 2 I s · A α k l s , T λ s I s , α k γ s I s , α k ( κ s + ω 0 s c s ) I s , α β p β 0 A 5 , k l s B k s ( κ s + ω 0 s c s ) p k Λ s , D a n t i s = 0 A 1 s A 2 s A 3 s A 4 s A 1 s 0 0 0 0 A 2 s 0 0 0 0 A 3 s 0 0 0 0 A 4 s 0 0 0 0 ,
where
A α β k s = ( p s , α δ β k + p s , β δ α k ) 2 3 δ α β p s , l δ l k , A α k l s , T = 2 ( p s , β δ α k δ β l p s , β 3 δ α β δ k l ) A 0 s = 2 η s I s , α k I s , β l + ( η ¯ s 2 3 η s ) δ k l I s , α β + α 1 s ( I s , α j Q s , j k I s , β l + I s , α k Q s , j l I s , β j ) + α 2 s Q s , k l I s , α i Q s , i j I s , β j + α 3 s ( I s , α i p s , i p s , k I s , β l + I s , α k p s , l p s , j I s , β j ) + α 4 s p s , k p s , l I s , α i p s , i p s , j I s , β j , B α s = λ 1 s p s , α ω p s · s p s , α ω p s , α s · p s , A 1 s = ν 0 s I s , α k I s , β l + a s [ I s , α j Q s , j k I s , β l + I s , α k Q s , j l I s , β j ] a s ( Q s , k l I s , α i ( Q s , i j + 1 3 δ i j ) I s , j β ) + θ 1 s δ k l I s , α β , A 2 s = ν 1 s 2 ( p s , j I s , α k I s , j β + I s , α j p s , j I s , β k ) + θ 2 s p s , k I s , α β , A 2 s = ν 1 s p s , β δ α k δ β l + θ 2 s p s , α δ k l , A 3 s = ν 2 s 2 ( p s , j I s , α k I s , j β + I s , α j p s , j I s , β k ) + θ 3 s p s , k I s , α β , A 3 s = ν 2 s p s , β I s , α k δ β l + θ 3 s I s , α β p s , β δ k l , A 4 s = I s , α i ( ζ s p s , i p s , j + ζ 0 s Q s , i j + β 1 s ( i p s , j + j p s , i ) + β 2 s p s · s Q s , i j + β 3 s s 2 Q s , i j ( ζ 2 s + ζ 3 s p s , i p s , i ) δ i j ) I s , j β , A 4 s = A 4 s δ α k δ β l , A 5 s = ζ 1 s ( p s p s p s 2 3 I ) + ( ζ 2 s ω 2 s s · p s ω 2 s p s · s ) Q s .
Notice that structures of the coefficient matrices in the bulk and on the surface are similar. The resultant transport equations on the boundary are summarized in the following:
t c s + s · ( ( v + ω 0 s μ ^ s I s · p s ) c s ) = s · ( γ s I s · s μ s + χ 2 s I s · G s · p s χ 2 s 3 t r ( G s ) I s · p s + λ s I s · h s + κ s I s · p s μ ^ s + ν 2 s I s · D s · p s + θ 3 s D s , γ γ I s · p s ) + j c , t p s + s · ( ( v + ω s μ ^ s p s ) p s ) s · v p s + Ω s · p s = 1 γ 1 s h s + χ 1 s G s · p s χ 1 s 3 t r ( G s ) p s λ s s μ s + λ 1 s μ ^ s p s + ν 1 s D s · p s + θ 2 s D s , γ γ p s , Q s t + s · ( ( v + ω 2 s μ ^ s p s ) Q s ) s · v Q s + Ω s · Q s · I s I s · Q s · Ω s a s [ I s · Q s · D s · I s + I s · D s · Q s · I s ] = 1 γ 2 s G s + χ 1 s 2 ( p s h s + h s p s ) χ 2 s 2 ( p s s μ s + s μ s p s ) + ζ 1 s μ ^ s p s p s + ν 0 s I s · D s · I s a s ( Q s : D s ) I s · Q s · I s + ( θ 1 s D s , γ γ a s Q s : D s 3 ) I s + ( χ 2 s p s · s μ s 3 χ 1 s p s · h s 3 μ ^ s ζ 1 s p s 2 3 ) I + ζ 2 s μ ^ s Q s , t ρ s + s · ( ρ s v ) = 0 , t ( ρ s v ) = s · ( σ s ρ s vv ) + F s , σ s s = σ s ( g c s μ s ) I s 1 2 I s · ( p s h s h s p s ) ( I s · Q s · I s · G s I s · G s · I s · Q s ) , F s + [ [ ( σ σ e ) ] ] · n + [ [ f μ c ] ] n K c s c s ( n · f c ) K p s p s · ( n · f p ) K Q s Q s : ( n · f Q ) 2 κ σ s · n = 0 , [ [ μ j ] ] · n + j c μ s 0 , σ s s = I s · ( Π s I s + 2 η s D s + α 1 s ( Q s · D s + D s · Q s · I s ) + α 2 s ( Q s : D s ) Q s · I s + α 3 s ( p s p s · D s + D s · p s p s · I s ) + α 4 s ( p s p s : D s ) p s p s · I s a s ( Q s · G s · I s + G s · Q s · I s ) ν 0 s G s · I s + a s ( Q s : G s ) Q s · I s ν 1 s 2 ( p s h s · I s + h s p s · I s ) + ν 2 s 2 ( p s s μ s · I s + s μ s p s · I s ) + ζ s μ ^ s p s p s · I s + ζ 0 s μ ^ s Q s · I s + β 1 s μ ^ s ( s p s + s p s T ) · I s + β 2 s μ ^ s p s · s Q s · I s + β 3 s μ ^ s s 2 Q s · I s ) ,
where Π s = ( η ¯ s 2 3 η s ) D s , γ γ a s Q s , k l G s , k l 3 + θ 1 s G γ γ + θ 2 s p s , γ h s , γ θ 3 s p s , γ s , γ μ s + ( ζ 2 s + ζ 3 s p s , γ p s , γ ) μ ^ s is the surface-excess pressure.
Let’s examine the inequality [ [ μ j ] ] · n + j c μ s 0 . If we assume j c = 0 and [ [ j ] ] · n = 0 , the inequality is satisfied. If on the other hand we assign [ [ j ] ] · n = 0 and j c = 1 γ s μ s , the inequality is satisfied. If in addition, we assume the flux j s as zero, the transport equation for c s is given by
t c s + s · ( c s v ) = 1 γ s μ s .
This is an Allen–Cahn type transport equation for the excessive surface concentration.
We have derived the governing system of equations for the bulk as well as for the free surface. The surface transport equations serve as the dynamic boundary conditions for the bulk transport equations at the free interface. The passive liquid crystal limit is embedded in the model for active liquid crystals. One important point we show in this derivation is that the boundary transport equations together with the bulk transport equations can be derived together through the generalized Onsager principle to ensure consistency in the model. If we impose that the passive limit is dissipative, the corresponding submatrix in the Onsager relation must be nonnegative definite. The model we present can be extended by adding additional higher order terms. The mathematical structure behind the model is persistent, in which variational structure and energy dissipation structure for the passive component are clearly delineated. We next reduce the model to the isotropic viscous fluid model so that one can examine the consistency in a simple fluid limit.

4.3. Static Boundary Conditions

In a simple case, we consider the static boundary conditions on the free surface. If we set γ 1 s , γ 2 s , γ s as zero, we recover the static boundary conditions given below,
μ s = 0 , h s = 0 , G s = 0 .
If g depends on the internal variables rather than their spatial gradients, these translate into the following.
K c f i c n i + g c s = 0 , K p f i p n i + g p s = 0 , K Q f i Q n i + g Q s = 0 ,
where we assume that outside of the active liquid crystal region is an isotropic viscous fluid, thus the jumps [ [ · ] ] across the surface are the limits of the variables to the surface in the active liquid crystal domain. These are most likely Robin type boundary conditions.
The kinematic boundary condition is given by the following
( t + v · ) Φ = 0 .
The kinetic boundary condition is
s · σ s s + s g + 2 κ g n + F s = 0 , F s = [ [ σ σ e ] ] · n [ [ f μ c ] ] n s c s g c s s p s · g p s s Q s : g Q s ,
where σ s s is defined in the above subsection, and we have used the boundary conditions (78) and σ s · n = 0 .
Next, we consider the specific apolar active liquid crystal system, whose model is given in Section 3. The bulk and surface free energy of the system are given by
f = f ( c , c , Q , Q ) , g = g 0 ( c s ) + g 1 Q s : nn .
Then
μ s = g 0 c s + K c f ( c ) · n , G s = g 1 ( nn ) + K Q f ( Q ) · n .
Thus the boundary conditions for c s , Q s on the free surface are μ s = 0 , G s = 0 . The kinetic boundary condition is
s · σ s s + s g + 2 κ g n + F s = 0 , F s = [ [ σ σ e ] ] · n [ [ f μ c ] ] n s c s g c s s Q s : g Q s , σ s s = I s · ( Π s I s + 2 η s D s + α 1 s ( Q s · D s + D s · Q s · I s ) + α 2 s ( Q s : D s ) Q s · I s + ζ 0 s μ ^ s Q s · I s + β 3 s μ ^ s s 2 Q s · I s ) .
The model we derived include many liquid crystal models and the viscous fluid model. For an isotropic viscous fluid, the surface stress tensor and elastic force are given by,
σ s = g I s , F s = [ [ σ ] ] · n .
The surface excessive momentum transport equation is given by
t ( ρ s v ) = s · ( σ s ρ s vv ) [ [ σ ] ] · n .
In the inertialess limit, the kinetic boundary condition reduces to
s · σ s [ [ σ ] ] · n = 0 .
This is the kinetic boundary condition commonly used in studies of fluid dynamics of viscous fluids with free surface boundaries.

5. Conclusions

We have derived systematically the transport equations together with the boundary conditions for an active liquid crystal in a confined fixed domain and in free surface domains using the generalized Onsager principle. The model reduces to a general hydrodynamic model for passive liquid crystals when the activity is absent, where energy dissipation is guarranteed. Likewise, the boundary conditions reduce to static boundary conditions. The resulting model equations incorporate many momentum-conserving models for liquid crystals with which to study complex interfacial phenomena involving active as well as passive liquid crystals. Detailed analyses of the model together with their applications to active matter systems and development of numerical schemes will be pursued in subsequent studies.

Acknowledgments

Jun Li’s work is partially supported by NSF of China through a grant (NSFC-11301287). M. Gregory Forest’s research is partially supported by NSF DMS-1517274 , AFOSR FA9550-12-1-0178, ARO-12-60317-MS. Qi Wang’s research is partially supported by the National Science Foundation through grants DMS-1200487 and DMS-1517347, AFOSR grant FA9550-12-1-0178, NIH grant NIH-R01GM078994-05A1, and an SC EPSCOR GEAR award.

Author Contributions

Xiaogang Yang, Jun Li, M. Gregory Forest and Qi Wang conceived and developed the theory; Xiaogang Yang and Jun Li performed the calculations and analyses; M. Gregory Forest and Qi Wang wrote the paper. All authors have read and approved the final manuscript.

Conflicts of Interest

The authors declare no conflict of interest.

Appendix A: Translational and Rotational Invariance of the Free Energy

In this appendix, we identify the Ericksen stress σ e from the translational invariance of the free energy and the antisymmetric stress σ a from the rotational invariance of the free energy. In addition, we use the angular momentum conservational law to prove symmetry of stress σ s .

A.1. Ericksen Stress σ e

Our goal is to drive the Ericksen stress σ e of this system. The free energy density f = f ( c , α c , p α , β p α , Q α β , γ Q α β ) . The total free energy is F = V d x f . If the volume V (with surface S = V ) of the system is changed by an amount δ V to volume V + δ V , the change of the free energy is
δ F = V + δ V d x f ( c + δ c , α c + δ α c , p α + δ p α , β p α + δ β p α , Q α β + δ Q α β , γ Q α β + δ γ Q α β ) V d x f ( c , α c , p α , β p α , Q α β , γ Q α β ) = δ V d x f ( c , α c , p α , β p α , Q α β , γ Q α β ) + V d x [ f c δ c + f ( α c ) δ α c + f p α δ p α + f ( β p α ) δ β p α + f Q α β δ Q α β + f ( γ Q α β ) δ γ Q α β ] = S d S α [ f u α ] + V d x [ μ δ c h α δ p α G α β δ Q α β ] + V d x [ α ( f ( α c ) δ c ) + β ( f ( β p α ) δ p α ) + γ ( f ( γ Q α β ) δ Q α β ) ] ( d x = u α d S α in δ V ) = S d S α [ f u α ] + V d x [ μ δ c h α δ p α G α β δ Q α β ] + S d S α [ f ( α c ) δ c ] + S d S β [ f ( β p α ) δ p α ] + S d S γ [ f ( γ Q α β ) δ Q α β ] = S d S α [ ( f μ c ) u α ] + V d x [ c α ( μ u α ) + h α ( β p α ) u β + G α β ( γ Q α β ) u γ ] S d S α [ f ( α c ) ( β c ) u β ] S d S β [ f ( β p α ) ( γ p α ) u γ ] S d S γ [ f ( γ Q α β ) ( θ Q α β ) u θ ] ,
where μ = δ F δ c = f c α ( f ( α c ) ) is the chemical potential, h α = δ F δ p α = f p α + β ( f ( β p α ) ) is the molecule field, G α β = δ F δ G α β = f Q α β + γ ( f ( γ Q α β ) ) is the conjugation field of Q .
Here, we have made use of δ c = ( α c ) u α and δ p α = ( β p α ) u β and δ Q α β = ( γ Q α β ) u γ , that means the movement of the matter volume is linear and there is no rotation, where u ( r ) is the change from V to V + δ V in the direction normal to the original boundary V . We obtain the expression of the Ericksen stress tensor σ e :
σ α β e = ( f c μ ) δ α β f ( β c ) α c f ( β p γ ) α p γ f ( β Q γ θ ) α Q γ θ .
In order to derive the corresponding Gibbs–Duhem relation, we first calculate the total differential of the free energy
d f = f c d c + f ( α c ) d ( α c ) + f p α d p α + f ( β p α ) d ( β p α ) + f Q α β d p α β + f ( β Q α γ ) d ( β Q α γ ) = μ d c h α d p α G α β d Q α β + α ( f ( α c ) d c ) + β ( f ( β p α ) d p α ) + γ ( f ( γ Q α β ) d Q α β ) .
Then, we have
d σ α β e = ( d f d ( c μ ) ) δ α β d ( f ( β c ) α c + f ( β p γ ) α p γ + f ( β Q γ θ ) α Q γ θ ) = [ c d μ h γ d p γ G γ θ d Q γ θ ] δ α β + γ ( f ( γ c ) d c ) δ α β + γ ( f ( γ p δ ) d p δ ) δ α β + γ ( f ( γ Q θ δ ) d Q θ δ ) δ α β d ( f ( β c ) α c + f ( β p γ ) α p γ + f ( β Q γ θ ) α Q γ θ ) .
It follows that
β σ α β e = [ c β μ h γ β p γ G γ θ β Q γ θ ] δ α β + γ ( f ( γ c ) β c ) δ α β + γ ( f ( γ p δ ) β p δ ) δ α β + γ ( f ( γ Q θ δ ) β Q θ δ ) δ α β β ( f ( β c ) α c + f ( β p γ ) α p γ + f ( β Q γ θ ) α Q γ θ ) = c α μ h γ α p γ G γ θ α Q γ θ .
This is the Gibbs–Duhem relation. Finally, the change of the free energy is given by
δ F = V d x [ c α ( μ u α ) + h α ( β p α ) u β + G α β ( γ Q α β ) u γ ] + S d S β [ σ α β e u α ] = V d x [ c α ( μ u α ) + h α ( β p α ) u β + G α β ( γ Q α β ) u γ + β ( σ α β e u α ) ] .
If the displacement u ( r ) is constant, that means the change is a pure linear translation, β u α = 0 , then
δ F = V d x [ c α μ + h γ ( α p γ ) + G β γ ( α Q β γ ) + β σ α β e ] u α = 0 .
This is the translational invariance of the free energy.

A.2. Antisymmetric Stress σ a

The free energy is also invariant under a pure uniform rotation around the original point with respect to an angle θ. The displacement u = θ × r , that is u α = ϵ α β γ θ β r γ . The changes of the variables are δ c = u β β c , δ p α = ( u β β p α + ϵ α β γ θ β p γ ) , δ Q α β = ( u γ γ Q α β + ϵ α δ γ θ δ Q γ β + ϵ β δ γ θ δ Q γ α ) . Under the pure uniform rotation, the convected co-rotational derivatives of a scale, a vector and a tensor are all zero. For the scalar, that is t c + v · c = 0 , then
δ c = c t δ t = v β β c δ t = u β β c ,
where the velocity is defined as v β = u β / δ t , that means that after a period of time δ t the displacement is u β . For a vector, we have t p + v · p + Ω · p = 0 , then
δ p α = p α t δ t = ( v β β p α + Ω α β p β ) δ t = u β β p α + ϵ α β γ θ β p γ ,
where Ω α β δ t = ϵ α γ β θ γ . For a tensor, the convected co-rotational derivative t Q + v · Q + Ω · Q Q · Ω = 0 , then
δ Q α β = Q α β t δ t = ( v γ γ Q α β + Ω α γ Q γ β Q α γ Ω γ β ) δ t = u γ γ Q α β + ϵ α δ γ θ δ Q γ β + ϵ β δ γ θ δ Q γ α .
We list some useful formulas as follows:
ϵ α β γ = ϵ γ β α , ϵ α α γ = 0 , β ( u β ) = ϵ α β γ θ β β ( r γ ) = ϵ α β γ θ β δ β γ = ϵ α β β θ β = 0 , δ ( u α ) = ϵ α β γ θ β δ ( r γ ) = ϵ α β γ θ β δ δ γ = ϵ α β δ θ β .
Then, we calculate the change of the free energy,
d F = S d S α [ f u α ] + V d x [ μ δ c h α δ p α G α β δ Q α β ] + V d x [ α ( f ( α c ) δ c ) + β ( f ( β p α ) δ p α ) + γ ( f ( γ Q α β ) δ Q α β ) ] = S d S α [ f c μ ] u α V d x [ β ( σ α β e ) u α ] V d x [ h α ϵ α β γ θ β p γ + G α β ( ϵ α δ γ θ δ Q γ β + ϵ β δ γ θ δ Q γ α ) ] + S d S β [ f ( β c ) α c f ( β p γ ) α p γ f ( β Q γ θ ) α Q γ θ ] u α + S d S δ [ f ( δ p α ) ϵ α β γ θ β p γ ] + S d S δ [ f ( δ Q α β ) ( ϵ α δ γ θ δ Q γ β + ϵ β δ γ θ δ Q γ α ) ] = S d S β [ σ α β e u α ] V d x [ β ( σ α β e ) u α ] V d x [ h α ϵ α β γ θ β p γ + G α β ( ϵ α δ γ θ δ Q γ β + ϵ β δ γ θ δ Q γ α ) ] + S d S δ [ f ( δ p α ) ϵ α β γ θ β p γ ] + S d S θ [ f ( θ Q α β ) ( ϵ α δ γ θ δ Q γ β + ϵ β δ γ θ δ Q γ α ) ] = V d x [ σ α β e ϵ α γ β θ γ ] V d x [ 1 2 ( h α p γ p α h γ ) + ( G α δ Q δ γ Q α δ G δ γ ) ] ϵ α β γ θ β + S d S δ [ f ( δ p α ) p γ + f ( δ Q α θ ) Q γ θ + f ( δ Q θ α ) Q γ θ ) ] ϵ α β γ θ β = V ϵ α β γ θ β [ σ α γ e + σ α γ a ] d x + S ϵ α β γ θ β p γ ( σ p · d S ) α + 2 S ϵ α β γ θ β [ ( Q · M Q · d S ) γ α ] = 0 ,
where we define σ α β a = 1 2 ( p α h β h α p β ) + ( Q α δ G δ γ G α δ Q δ γ ) , σ α β p = f ( β p α ) , M α β γ Q = f ( γ Q α β ) . M α β γ Q is symmetric about α , β . Because the rotation angle is arbitrary, we get the relationship of the angular momentum.
V ϵ α β γ [ σ α γ e + σ α γ a ] d x + S ϵ α β γ p γ ( σ p · d S ) α + 2 S ϵ α β γ [ ( Q · M Q · d S ) γ α ] = 0 .
We show the z-component clearly,
V d x { [ σ y x e σ x y e ] + [ σ y x a σ x y a ] } + S [ p × ( σ p · d S ) ] z + 2 S ϵ γ z α [ ( Q · M Q · d S ) γ α ] = 0 .
The second term on the right-hand side is the torque on the surface due to polar order p , and the last term is the torque on the surface due to nematic order tensor Q .

A.3. Angular Momentum Conservation Law

We have calculated that the dissipation of the total energy
d E t o t o l d t = V d x ( β v α ) σ α β s + j α ( α μ ¯ ) + P α h α + G α β Q α β + r μ ^ ,
where σ α β s = σ α β σ α β e σ α β a . We need to prove that σ s is symmetry using the angular momentum conservation. The angular momentum conservation law is given by
t V d x ( r × ρ v ) = S r × ( ( σ ρ vv ) · d S ) + S p × ( σ p · d S ) + 2 S ϵ α β γ [ ( Q · M Q · d S ) γ α ] ,
where σ α β p = f ( β p α ) and M α β γ Q = f ( γ Q α β ) are defined above. The term on the left-hand side is the derivative of the total angular momentum in the volume. The terms on the right-hand side are the total torque acting on the surface. There are three contributions: the first one is the torque of the stress ( σ ρ vv ) acting on the surface; the second one is the torque of the polarity vector ( p ) acting on the surface and the last one is the torque of the nematic order tensor ( Q ) acting on the surface. We use the linear momentum conservation law to transform the left-hand side of Equation (A16) and use the Equation (A14) to arrive at:
V d x [ σ x y σ y x ] = V d x { [ σ y x e σ x y e ] + [ σ y x a σ x y a ] } ,
where we only show the z-component. Analogously, we can get the x , y -components. Because σ a is antisymmetric, then we get the antisymmetric part
1 2 [ ( σ σ e ) α β ( σ σ e ) β α ] = σ α β a .
So, the stress σ α β s = σ α β σ α β a σ α β e is symmetric, i.e., β v α σ α β s = D : σ s .

A.4. The Convected Co-Rotational Derivatives:

The convected co-rotational derivatives of the polarization and nematic order are defined by
p ˙ = p t + v · p + Ω · p , Q = Q t + v · Q + [ Ω · Q Q · Ω ] ,
where v is the velocity of the background fluid, Ω α β = 1 2 ( α v β β v α ) is the vorticity tensor. The vector p ˙ and tensor Q represent the rate of change of the field p and Q with respect to the background fluid, which are translational and rotational invariant.

A.4.1. Translational Invariance

If both the background fluid and the director p move with a constant velocity v = v 0 , then Ω = 0 , Ω · p = 0 . At time t, the director at the location r is p ( r , t ) . After a small time Δ t , at time t + Δ t , the director at the location r + v 0 Δ t is p ( r + v 0 Δ t , t + Δ t ) = p ( r , t ) .
p t = lim Δ t 0 p ( r , t + Δ t ) p ( r , t ) Δ t = lim Δ t 0 p ( r , t + Δ t ) p ( r + v 0 Δ t , t + Δ t ) + p ( r + v 0 Δ t , t + Δ t ) p ( r , t ) Δ t = lim Δ t 0 v 0 · p ( r + a v 0 Δ t , t + Δ t ) Δ t Δ t = v · p ,
where a is some constant, 0 < a < 1 . So, p ˙ = p t + v · p + Ω · p = 0 under a pure linear translation. Similarly, we can show Q = Q t + v · Q + [ Ω · Q Q · Ω ] = 0 .

A.4.2. Rotational Invariance

If both the background fluid and the director p rotate around a certain axis at the same angular velocity θ, velocity v = θ × r , and Ω · p = θ × p . At time t, the director at the location r is p ( r , t ) . After a small time Δ t , at time t + Δ t , the director at the location r + θ × r Δ t is p ( r + θ × r Δ t , t + Δ t ) = p ( r , t ) + θ × p ( r , t ) Δ t .
p t = lim Δ t 0 p ( r , t + Δ t ) p ( r , t ) Δ t = lim Δ t 0 p ( r , t + Δ t ) p ( r + θ × r Δ t , t + Δ t ) + p ( r + θ × r Δ t , t + Δ t ) p ( r , t ) Δ t = lim Δ t 0 ( θ × r ) · p ( r + ( θ × r ) a Δ t , t + Δ t ) Δ t + θ × p ( r , t ) Δ t Δ t = ( θ × r ) · p ( r , t ) + θ × p ( r , t ) = v · p Ω · p ,
where a is some constant, 0 < a < 1 . So, p ˙ = p t + v · p + Ω · p = 0 under a pure rotation. For tensor Q , we recall from its definition Q = V ( ν ^ ν ^ I d ) f d ( r , ν ^ , t ) d ν ^ = M c I d . Using the result of p ,
M ( r + θ × r Δ t , t + Δ t ) = V ( ν ^ + θ × ν ^ ) ( ν ^ + θ × ν ^ Δ t ) f d ( r , ν ^ , t ) d ν ^ = M ( r , t ) + V ( ν ^ θ × ν ^ + θ × ν ^ ν ^ ) f d ( r , ν ^ , t ) d ν ^ Δ t + O ( Δ 2 t ) = M ( r , t ) V ( ν ^ [ Ω · ν ^ ] + [ Ω · ν ^ ] ν ^ ) f d ( r , ν ^ , t ) d ν ^ Δ t + O ( Δ 2 t ) .
The time derivative is
M t = lim Δ t 0 M ( r , t + Δ t ) M ( r , t ) Δ t = lim Δ t 0 M ( r , t + Δ t ) M ( r + θ × r Δ t , t + Δ t ) + M ( r + θ × r Δ t , t + Δ t ) M ( r , t ) Δ t = lim Δ t 0 ( θ × r ) · M ( r + ( θ × r ) a Δ t , t + Δ t ) Δ t V ( ν ^ [ Ω · ν ^ ] + [ Ω · ν ^ ] ν ^ ) f d ( r , ν ^ , t ) d ν ^ Δ t + O ( Δ 2 t ) Δ t = ( θ × r ) · M ( r , t ) V ( ν ^ [ Ω · ν ^ ] + [ Ω · ν ^ ] ν ^ ) f d ( r , ν ^ , t ) d ν ^ = v · M + M · Ω Ω · M .
Notice Q = M c I d , we have the rotational invariance on Q : Q = Q t + v · Q + [ Ω · Q Q · Ω ] = 0 .

Appendix B: Surface Divergence Theorem

The surface divergence theorem is [53]:
Γ d l i R j k l , . . . = Γ [ s , i R j k l , . . . + 2 κ n i R j k l , . . . ] d s = Γ s , α · ( I s , α i R j k l , . . . ) d s .
This leads to the integration by parts formula,
Γ [ s , i ( R j k l , . . . ) S α β , γ , . . . ] d s = Γ R j k l , . . . s , i S α β , γ , . . . d s + Γ d l i R j k l , . . . S α β , γ , . . . Γ 2 κ n i R j k l , . . . S α β , γ , . . . d s ,
where κ = 1 2 ( κ 1 + κ 2 ) is the mean surface curvature, κ 1 , κ 2 are the principal curvatures of the surface.

References

  1. Marchetti, M.; Joanny, J.; Ramaswamy, S.; Liverpool, T.; Prost, J.; Rao, M.; Simha, R.A. Hydrodynamics of soft active matter. Rev. Mod. Phys. 2013, 85. [Google Scholar] [CrossRef]
  2. Ramaswamy, S. The Mechanics and Statistics of Active Matter. Annu. Rev. Condens. Matter Phys. 2010, 1, 323–345. [Google Scholar] [CrossRef]
  3. Baskaran, A.; Marchetti, M. Self-regulation in self-propelled nematic fluids. Eur. Phys. J. E 2012, 35, 1–8. [Google Scholar] [CrossRef] [PubMed]
  4. Aditi, S.R.; Ramaswamy, S. Hydrodynamic Fluctuations and Instabilities in Ordered Suspensions of Self-Propelled Particles. Phys. Rev. Lett. 2002, 89, 058101. [Google Scholar] [CrossRef] [PubMed]
  5. De Gennes, P.G.; Prost, J. The Physics of Liquid Crystals; Oxford Science Publications: Oxford, UK, 1993. [Google Scholar]
  6. Liverpool, T.B.; Marchetti, M.C. Rheology of Active Filament Solutions. Phys. Rev. Lett. 2006, 97. [Google Scholar] [CrossRef] [PubMed]
  7. Gruler, H.; Dewald, U.; Eberhardt, M. Nematic liquid crystals formed by living amoeboid cells. Eur. Phys. J. B 1999, 11, 187–192. [Google Scholar] [CrossRef]
  8. Kemkemer, R.; Kling, D.; Kaufmann, D.; Gruler, H. Elastic properties of nematoid arrangements formed by amoeboid cells. Eur. Phys. J. E 2000, 1, 215–225. [Google Scholar] [CrossRef]
  9. Saintillan, D.; Shelley, M.J. Instabilities and Pattern Formation in Active Particle Suspensions: Kinetic Theory and Continuum Simulations. Phys. Rev. Lett. 2008, 100. [Google Scholar] [CrossRef] [PubMed]
  10. Forest, M.G.; Phuworawong, P.; Wang, Q.; Zhou, R. Rheological signatures in limit cycle behavior of dilute, active, polar LCPs in steady shear. Philos. Trans. R. Soc. A 2014. [Google Scholar] [CrossRef] [PubMed]
  11. Forest, M.G.; Zhou, R.; Wang, Q. Kinetic attractor phase diagrams of active nematic suspensions: The dilute regime. Soft Matter 2015, 11, 6393–6402. [Google Scholar] [CrossRef] [PubMed]
  12. Vicsek, T.; Czirók, A.; Ben-Jacob, E.; Cohen, I.; Shochet, O. Novel Type of Phase Transition in a System of Self-Driven Particles. Phys. Rev. Lett. 1995, 75, 1226–1229. [Google Scholar] [CrossRef] [PubMed]
  13. Baskaran, A.; Marchetti, M. Enhanced Diffusion and Ordering of Self-Propelled Rods. Phys. Rev. Lett. 2008, 101. [Google Scholar] [CrossRef] [PubMed]
  14. Baskaran, A.; Marchetti, M. Statistical mechanics and hydrodynamics of bacterial suspensions. Proc. Natl. Acad. Sci. USA 2009, 106, 15567–15572. [Google Scholar] [CrossRef] [PubMed]
  15. Grégoire, G.; Chaté, H. Onset of Collective and Cohesive Motion. Phys. Rev. Lett. 2004, 92. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  16. Chaté, H.; Ginelli, F.; Grégoire, G.; Raynaud, F. Collective motion of self-propelled particles interacting without cohesion. Phys. Rev. E 2008, 77. [Google Scholar] [CrossRef] [PubMed]
  17. Ginelli, F.; Peruani, F.; Bär, M.; Chaté, H. Large-Scale Collective Properties of Self-Propelled Rods. Phys. Rev. Lett. 2010, 104. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  18. Peshkov, A.; Aranson, I.S.; Bertin, E.; Chaté, H.; Ginelli, F. Nonlinear Field Equations for Aligning Self-Propelled Rods. Phys. Rev. Lett. 2012, 109. [Google Scholar] [CrossRef] [PubMed]
  19. Nagai, K.H.; Sumino, Y.; Montagne, R.; Aranson, I.S.; Chaté, H. Collective Motion of Self-Propelled Particles with Memory. Phys. Rev. Lett. 2015, 114. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  20. Bertin, E.; Droz, M.; Grégoire, G. Boltzmann and hydrodynamic description for self-propelled particles. Phys. Rev. E 2006, 74, 022101. [Google Scholar] [CrossRef] [PubMed]
  21. Bertin, E.; Droz, M.; Grégoire, G. Hydrodynamic equations for self-propelled particles: Microscopic derivation and stability analysis. J. Phys. A: Math. Theor. 2009, 42. [Google Scholar] [CrossRef]
  22. Bertin, E.; Chaté, H.; Ginelli, F.; Mishra, S.; Peshkov, A.; Ramaswamy, S. Mesoscopic theory for fluctuating active nematics. New J. Phys. 2013, 15, 085032. [Google Scholar] [CrossRef]
  23. Mishra, S.; Baskaran, A.; Marchetti, M.C. Fluctuations and pattern formation in self-propelled particles. Phys. Rev. E 2010, 81, 061916. [Google Scholar] [CrossRef] [PubMed]
  24. Saintillan, D.; Shelley, M.J. Orientational Order and Instabilities in Suspensions of Self-Locomoting Rods. Phys. Rev. Lett. 2007, 99, 058102. [Google Scholar] [CrossRef] [PubMed]
  25. Saintillan, D.; Shelley, M.J. Instabilities, pattern formation, and mixing in active suspensions. Phys. Fluids 2008, 20. [Google Scholar] [CrossRef]
  26. Hohenegger, C.; Shelley, M.J. Stability of active suspensions. Phys. Rev. E 2010, 81, 046311. [Google Scholar] [CrossRef] [PubMed]
  27. Baskaran, A.; Marchetti, M.C. Nonequilibrium statistical mechanics of self-propelled hard rods. J. Stat. Mech. 2010, 2010, 04019. [Google Scholar] [CrossRef]
  28. Liverpool, T.B.; Marchetti, M.C.; Joanny, J.F.; Prost, J. Mechanical response of active gels. Europhys. Lett. 2009, 85. [Google Scholar] [CrossRef]
  29. Gopinath, A.; Hagan, M.; Marchetti, M.; Baskaran, A. Dynamical self-regulation in self-propelled particle flows. Phys. Rev. E 2012, 85, 061903. [Google Scholar] [CrossRef] [PubMed]
  30. Forest, M.G.; Zhou, R.; Wang, Q. A kinetic theory and its predictions for semidilute active nematic suspensions. Soft Matter 2013, 9, 5207–5222. [Google Scholar] [CrossRef]
  31. Baskaran, A.; Marchetti, M. Hydrodynamics of self-propelled hard rods. Phys. Rev. E 2008, 77, 011920. [Google Scholar] [CrossRef] [PubMed]
  32. Wang, Q.; Yang, X.; David, A.; Elston, T.; Jacobson, K.; Maria, M.; Forest, M.G. Computational and Modeling Strategies for Cell Motility: Bacteria as Multicellular Organisms; Springer: New York, NY, USA, 2012; pp. 257–296. [Google Scholar]
  33. Kruse, K.; Joanny, J.; Jülicher, F.; Prost, J.; Sekimoto, K. Asters, Vortices and Rotating Spirals in Active Gels of Polar Filaments. Phys. Rev. Lett. 2004, 92, 078101. [Google Scholar] [CrossRef] [PubMed]
  34. Voituriez, R.; Joanny, J.; Prost, J. Spontaneous flow transition in active polar gels. Europhys. Lett. 2005, 70. [Google Scholar] [CrossRef]
  35. Joanny, J.; Jülicher, F.; Kruse, K.; Prost, J. Hydrodynamic theory for multi-component active polar gels. New J. Phys. 2007, 9. [Google Scholar] [CrossRef]
  36. Yang, X. Modeling and Numerical Simulations of Active Liquid Crystals. Ph.D. Thesis, Nankai University, Tianjin, China, 2014. [Google Scholar]
  37. Onsager, L. Reciprocal Relations in Irreversible Processes. I. Phys. Rev. 1931, 37, 405–426. [Google Scholar] [CrossRef]
  38. Onsager, L. Reciprocal Relations in Irreversible Processes. II. Phys. Rev. 1931, 38, 2265–2279. [Google Scholar] [CrossRef]
  39. Zhang, J.; Xu, X.; Qian, T. Anisotropic particle in viscous shear flow: Navier slip, reciprocal symmetry, and Jeffery orbit. Phys. Rev. E 2015, 91, 033016. [Google Scholar] [CrossRef] [PubMed]
  40. De Groot, S.; Mazur, P. Non-Equilibrium Thermodynamics; Dover: New York, NY, USA, 1984. [Google Scholar]
  41. Kondepudi, D.; Prigogine, l. Modern Thermodynamics From Heat Engines to Dissipative Structures; John Wiley & Sons: New York, NY, USA, 1998. [Google Scholar]
  42. Lebon, G.; Jou, D.; Casas-Vazquez, J. Understanding Non-equilibrium Thermodynamics–Foundations, Applications, Frontiers; Springer-Verlag: Berlin/Heidelberg, Germany, 2008. [Google Scholar]
  43. Müller, I. A History of Thermodynamics The Doctrine of Energy and Entropy; Springer-Verlag: Berlin/Heidelberg, Germany, 2007. [Google Scholar]
  44. Beris, A.N.; Edwards, B. Thermodynamics of Flowing Systems with Internal Microstructure; Oxford University Press: Oxford, UK, 1994. [Google Scholar]
  45. Ottinger, H.C. Beyond Equilibrium Thermodynamics; Wiley Online Library: New York, NY, USA, 2006. [Google Scholar]
  46. Sagis, L.M. Generic model for multiphase systems. Adv. Colloid Interface Sci. 2010, 153, C58–C69. [Google Scholar] [CrossRef] [PubMed]
  47. Sagis, L.M. Dynamic properties of interfaces in soft matter: Experiments and theory. Rev. Modern Phys. 2011, 83, 1367–1403. [Google Scholar] [CrossRef]
  48. Sagis, L.M. Rheology of interfaces stabilized by a 2D suspension of anisotropic particles: A classical irreversible thermodynamics theory. Soft Matter 2011, 7, 7727–7736. [Google Scholar] [CrossRef]
  49. Sagis, L.M. Dynamic behavior of interfaces: Modeling with nonequilibrium thermodynamics. Adv. Colloid Interface Sci. 2014, 206, 328–343. [Google Scholar] [CrossRef] [PubMed]
  50. Jou, D.; Casas-Vhzquez, J.; Lebon, G. Extended irreversible thermodynamics. Rep. Prog. Phys. 1988, 51, 1105–1179. [Google Scholar] [CrossRef]
  51. Chandrasekhar, S. Liquid Crystals; Cambridge University Press: Cambridge, UK, 1992. [Google Scholar]
  52. Leslie, F.M. Theory of Flow Phenomena in Liquid Crystals; Academic Press: New York, NY, USA, 1979. [Google Scholar]
  53. Edwards, D.; Brenner, H.; Wasan, D. Interfacial Transport Processes and Rheology; Butterworth-Heinemann: Oxford, UK, 1991. [Google Scholar]
  54. Slattery, J.C.; Sagis, L.; Oh, E.-S. Interfacial Transport Phenomena, 2nd ed.; Springer: New York, NY, USA, 2007. [Google Scholar]

Share and Cite

MDPI and ACS Style

Yang, X.; Li, J.; Forest, M.G.; Wang, Q. Hydrodynamic Theories for Flows of Active Liquid Crystals and the Generalized Onsager Principle. Entropy 2016, 18, 202. https://doi.org/10.3390/e18060202

AMA Style

Yang X, Li J, Forest MG, Wang Q. Hydrodynamic Theories for Flows of Active Liquid Crystals and the Generalized Onsager Principle. Entropy. 2016; 18(6):202. https://doi.org/10.3390/e18060202

Chicago/Turabian Style

Yang, Xiaogang, Jun Li, M. Gregory Forest, and Qi Wang. 2016. "Hydrodynamic Theories for Flows of Active Liquid Crystals and the Generalized Onsager Principle" Entropy 18, no. 6: 202. https://doi.org/10.3390/e18060202

APA Style

Yang, X., Li, J., Forest, M. G., & Wang, Q. (2016). Hydrodynamic Theories for Flows of Active Liquid Crystals and the Generalized Onsager Principle. Entropy, 18(6), 202. https://doi.org/10.3390/e18060202

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop