Mathematics and liquid crystals

ABSTRACT

A review is given of some mathematical contributions, ideas and questions concerning liquid crystals.

1. Introduction

Recent years have seen a remarkable increase in the number of mathematicians¹ working on problems of liquid crystals. This reflects both the increased ability of modern mathematics to provide useful information on the problems of liquid crystals, and the interest and challenges of these problems for mathematicians. Indeed the questions raised by liquid crystals are generating new mathematical techniques, which can be expected in turn to influence applications of mathematics to other parts of science.

Everyone who works on liquid crystals uses mathematics in one form or another, and professional mathematicians (for example, those who work in mathematics departments) don’t have a prerogative on how mathematics should best be used. Nevertheless mathematicians can have a different way of viewing things that adds new perspectives, and it is one aim of this article to convey something of this. The ultimate objective should be for mathematics to illuminate experiments and phenomena, to help predict and interpret experimental results, and to help suggest new experiments. Of course the hypotheses of mathematical results may often ignore relevant physical effects. However by making simplifying assumptions it may be possible to prove a rigorous theorem that both gives insight into what more complex models may predict, resolving ambiguities arising from other theoretical predictions based on approximations or numerical computation, and enabling placing of the problem within the existing mathematical body of knowledge.

This article is intended less for mathematicians than for other researchers who may be interested in what mathematics can contribute to the field. Whereas I have attempted to give clear statements of results, at some points these are more informal than would be the case for a mathematical audience.

2. The Oseen-Frank model and the description of defects

Consider a nematic liquid crystal at rest at constant temperature in the absence of applied electromagnetic fields. (Although most of the applications of liquid crystals to technology depend on their interaction with electromagnetic fields, for simplicity such interactions are not discussed in this paper, since many of the key analytic difficulties already arise in the case when such fields are absent.)

Suppose that the liquid crystal completely fills a container , which we assume to be a bounded open set with boundary ∂Ω. The configuration of the liquid crystal at the point x is described by the director n(x), a unit vector giving the mean orientation at x of the rod-like molecules comprising the liquid crystal. Thus n(x) ∈ S², where S² denotes the the unit sphere of .

2.1. The Oseen-Frank energy

In the Oseen-Frank model the total free energy is given by (2.1) where (2.2) Here K₁, …, K₄ are the Frank constants, which are usually assumed to satisfy the strict form of the Ericksen inequalities [58] (2.3) which are necessary and sufficient that (2.4) for all n and some constant c₀ > 0. The first three terms in (2.2(2.2) ) describe respectively splay, twist and bend of the director field. The fourth saddle-splay term is a null Lagrangian, its integral depending only on the values of n on the boundary ∂Ω. Thus if n|_∂Ω is prescribed, such as for homeotropic or planar boundary conditions, the term can be ignored. However this is not the case if n|_∂Ω is only partially prescribed, such as for the planar degenerate boundary condition on part of the boundary having normal (see, for example, [106]), or for weak anchoring boundary conditions.

An important identity is (2.5) So if K₁ = K₂ = K₃, K₄ = 0 (the one-constant approximation) then (2.6) which is the energy functional for harmonic maps.

2.2. Orientability

The Oseen-Frank theory regards n as a vector field. However, due to the statistical head-to-tail symmetry of the constituent molecules, from a physical point of view n(x) is indistinguishable from − n(x). Hence n(x) is better thought of as a line field, or equivalently as a map from Ω to the set of all lines through the origin. A typical such line with direction ± n can be represented by the matrix n⊗n having components n_in_j. The set of such lines forms the real projective plane . Continuous (even smooth) line fields need not be orientable, so that it is impossible to assign a direction that turns them into a continuous vector field (see Fig. 1). We return later to the implications of this for the Oseen-Frank theory. For a detailed discussion see [21, 22].

Figure 1. A smooth line field in the exterior of a cylinder that cannot be oriented to make it a continuous vector field. The line field is parallel to the curves shown, with zero component in the perpendicular direction, so that the problem is two-dimensional. Indicated is an attempt to orient the line field, leading to a conflict in the ellipse shown.

2.3. Energy minimization and the Euler-Lagrange equation

The fundamental energy minimization problem is to find n that minimizes subject to the unit vector constraint |n| = 1 and suitable boundary conditions, for example , where is given. This is a problem of the calculus of variations.

If n is a minimizer and is any smooth mapping vanishing in a neighbourhood of ∂Ω then, for sufficiently small |ϵ|, (2.7) satisfies |n_ϵ(x)| = 1 and . Since I(n_ϵ) is minimized at ϵ = 0 we have that , provided this derivative exists. Noting that we obtain the weak form of the Euler-Lagrange equation, that for all such m or, in components, where repeated indices are summed from 1 to 3. Hence, integrating by parts and using the arbitrariness of m, we formally obtain the Euler-Lagrange equation a system of second order nonlinear PDE that is to be solved subject to the pointwise constraint |n| = 1.

This can be written in the equivalent form (2.8) where λ(x) is a Lagrange multiplier, or in components In the one-constant approximation (EL) becomes the harmonic map equation (2.9) How can we solve these equations? Are there any exact solutions?

2.4. Universal solutions

The question of what n(x) can be solutions of (EL) (smooth in some region of ) for all K₁, K₂, K₃, K₄, so called universal solutions, was addressed by Marris [130, 131], following Ericksen [59]. Marris showed that these consist of

(i)	constant vector fields, or those orthogonal to families of concentric spheres or cylinders,
(ii)	pure twists, such as (2.10)
(iii)	planar fields that form concentric or coaxial circles.

An example from family (i) is the hedgehog (2.11) which represents a point defect. Of course is not even continuous at 0, but for x ≠ 0 it is smooth and we have so that formally calculating its energy over the ball , noting that W(n, ∇n) ⩽ C|∇n|² for some constant C > 0, we find that More precisely, belongs to the Sobolev space (2.12) For a precise definition of Sobolev spaces the reader is referred to standard texts on partial differential equations, such as [63]. See also the discussion in [21]. The most important point in giving a precise definition is that ∇n has to be defined in a suitable weak sense. For the purposes of this article one can informally think of H¹(Ω; S²) as the space of finite-energy configurations for the Oseen-Frank theory; it is an example of a function space, that specifies allowed singularities in a mathematical model.

2.5. Existence of solutions and their regularity

A routine use of the so-called ‘direct method’ of the calculus of variations gives:

Theorem 1.

Let K₁ > 0, K₂ > 0, K₃ > 0, and let . Then there exists n* that minimizes I(n) over all n ∈ H¹(Ω; S²) with , and any such minimizer n* satisfies (WEL).

Remark 1.

Because n is specified on the whole of ∂Ω and because the saddle-splay term is a null Lagrangian, we do not need any conditions on K₄ to get the existence of a minimizer. However, for other boundary conditions such as planar degenerate we would need to assume the Ericksen inequalities (2.3(2.3) ).

Theorem 1 just tells us that there is some energy-minimizing director configuration n*, that is I(n) ⩾ I(n*) for all n with , and that any minimizer satisfies (WEL), and not what minimizers look like. Nevertheless the conclusion of the theorem is not obvious. Nor is it obvious that just because the problem posed has a sound physical basis the existence of a minimizer is necessarily assured, examples to the contrary arising, for example, in models of martensitic phase transformations [17], where nonattainment of minimum energy leads to an understanding of the appearance of very fine microstructures.

In the one-constant approximation (2.6(2.6) ) with K₁ > 0 there are deep results asserting more.

Theorem 2

(Schoen & Uhlenbeck [159], Brezis, Coron & Lieb [34]). In the one-constant approximation with K₁ > 0 any minimizer n* is smooth in Ω except for a finite number of point defects located at points x(i) ∈ Ω, and for some .

Theorem 3

(Brezis, Coron & Lieb [34]). In the one-constant approximation the hedgehog minimizes I subject to its own boundary conditions.

There is an elegant alternative proof of Theorem 3 due to Lin [118]. The discovery of Theorems 2 and 3 was motivated in particular by experiments concerning point defects and their annihilation [175], and by corresponding analytical and numerical studies [48, 80, 82, 83].

What about other solutions to (WEL) (these are often called weak solutions to (EL))? Unfortunately there seem to be far too many. In the one-constant approximation Rivière [153] showed that for a domain Ω with smooth boundary, with smooth on ∂Ω and not constant, there are infinitely many weak solutions, and that they can be discontinuous everywhere. The idea behind the construction is to insert dipoles (pairs of point defects with opposite topological charge) into a non constant smooth map. This is reminiscent of other physical systems of nonlinear partial differential equations where there are too many solutions and no satisfactory way of selecting the physical one, such as exhibited by Scheffer [157], Shnirelman [163] and De Lellis & Székelyhidi [55] for the Euler equations of inviscid flow.

The derivation of (WEL) used the outer variation (2.7(2.7) ). However, one can also take the inner variation (2.13) where is a smooth function vanishing in a neighbourhood of ∂Ω, and z_ϵ(x) is implicitly defined by (2.14) This variation rearranges the values of n in Ω and preserves the unit vector constraint. For sufficiently small |ϵ| (2.14(2.14) ) has a smooth solution z_ϵ: Ω → Ω with z_ϵ(x) = x near ∂Ω. Changing variables from x to z = z_ϵ we see that (2.15) Since I(n_ϵ) is minimized at ϵ = 0 we have that , giving (2.16) for all such . This is a weak form of the equation (2.17) If n is a smooth solution of (EL) then it is easily verified that n also satisfies (2.17(2.17) ). However solutions of (WEL) do not necessarily satisfy (2.16(2.16) ). Solutions of (WEL) that also satisfy (2.16(2.16) ) are called stationary. In the one-constant approximation it was proved by Evans [62] that any stationary solution is partially regular, that is such a solution is smooth outside a closed subset E of Ω having zero one-dimensional Hausdorff measure. Thus E cannot, for example, contain a line segment or arc. In view of the result of Rivière mentioned above, this means that stationary solutions of (WEL) are more regular than general solutions.

2.6. Dynamics and weak equilibrium solutions

It is interesting to examine the dynamical significance of the different notions of weak solution². The relevant dynamical equations are the Ericksen-Leslie equations (see [57, 115, 116] and for a modern treatment [164, Chapter 3]). For an incompressible fluid in the absence of body forces and couples, these equations for the velocity v(x, t) and director n(x, t) take the form (2.18) (2.19) (2.20) where ρ is the constant density, T is the Cauchy stress tensor, is the dissipation potential, , and . For the purposes of this discussion all we need to know is that when v = 0 and n = n(x), (2.21) where p is the pressure (a Lagrange multiplier corresponding to the incompressibility constraint (2.19(2.19) )). Note that it is indeed the dynamical equations for the flow of an incompressible liquid crystal that are appropriate here. Were we to consider compressible flow then the corresponding free-energy density W would depend not only on n, ∇n but also on the density ρ and perhaps ∇ρ (see [164, Chapter 5]), leading to a different variational problem for a functional I(ρ, n).

The requirement that v = 0, n = n(x), p = p(x) (where we assume that p is integrable over Ω) is a weak solution of (2.18(2.18) )–(2.20(2.20) ) is then that (WEL) holds and that (2.22) for all smooth vanishing in a neighbourhood of ∂Ω. If n is a stationary solution of (WEL) then (2.22(2.22) ) holds with p = p₀ − W, where p₀ is constant. However, it is not clear whether or not a solution n of (WEL) satisfying (2.22(2.22) ) is stationary. This would follow if we knew that n were partially regular, but this is also unclear.

It is interesting (see Ericksen [60]) that (2.18(2.18) ) can admit solutions with v = 0 but n = n(x, t), p = p(x, t), for which similar considerations apply.

2.7. Line defects

Another solution of type (i) is the two-dimensional hedgehog given by (2.23) This has the line defect . However, for any bounded domain Ω containing part of L, we have that some cylinder is contained in Ω, so that using (2.4(2.4) ) (2.24) Thus the two-dimensional hedgehog has infinite energy, and .

Other line defects are given by index defects, such as that illustrated in Fig. 2(a). Here the director lies in the (x₁, x₂) −plane and is parallel to the curves shown, with zero component in the x₃ direction. Such configurations cannot be described by the Oseen-Frank theory for two reasons. First, they are not orientable (see Fig. 2(b), where an attempt to orient it is illustrated, giving a conflict in the elliptical cylinder shown) so that they cannot be represented by a vector field n that is continuous outside the defect. Second, even if we restrict attention to one sector in which the line-field is orientable (such as the shaded region in Fig. 2(c)) the corresponding energy is infinite.

Figure 2. An index defect.

The problem that and the index singularities have infinite energy could potentially be fixed by modifying the growth of W(n, ∇n) for large |∇n| to be subquadratic, i.e. W(n, ∇n) ⩽ C(1 + |∇n|^p) for some C > 0, 1 ⩽ p < 2. After all, the quadratic dependence of W in ∇n suggests a theory designed to apply for small values of |∇n|, whereas near defects |∇n| is very large. As described in [16] it is easy to modify W to have subquadratic growth without affecting the Frank constants or other desirable properties. We return in Section 3 to the question of how the lack of orientability of the index defect might be handled.

3. The Lavrentiev phenomenon and function spaces

The Lavrentiev phenomenon, discovered in 1926 by Lavrentiev [113], is of profound importance for mechanics and physics. It occurs even in harmless looking problems of the one-dimensional calculus of variations³. It can be expressed as the statement:

3.1. An example from solid mechanics

For physical systems the Lavrentiev phenomenon has the uncomfortable implication that the function space is part of the model. An interesting illustration of this in a physical problem comes from nonlinear elasticity theory. Consider a nonlinear elastic body (composed of rubber, say) occupying in a reference configuration the unit ball . A typical deformation is described by a map , where y(x) denotes the deformed position of the material point x ∈ B. We look for a minimizer y of the total elastic free energy (3.1) subject to the boundary condition (3.2) where λ > 0, corresponding to a uniform radial expansion at the boundary (see Fig. 3(a)). We take as an example the compressible neo-Hookean material having free-energy density (3.3) where μ > 0 and h: (0, ∞) → [0, ∞) is a convex function satisfying (3.4) which is a standard model (though not the best) for rubber.

Figure 3. Solutions to the radial expansion problem for an elastic ball in different function spaces. (a) shows the boundary conditions, (b) the minimizer in the class of smooth maps, (c) the minimizer among radial maps in , with a cavity at the origin, (d) a possible minimizer in the space .

Among smooth mappings y the unique minimizer of I(y) subject to (3.2(3.2) ) is given (see Fig. 3(b)) by the uniform dilatation y(x) = λx (this follows from the fact that, due to the convexity of h, W is polyconvex, hence quasiconvex; for the explanation of these terms and more details see [14]). However, if λ is sufficiently large, the minimizer in the Sobolev space is not given by a uniform dilatation, because it is energetically more convenient to create a cavity (see Fig. 3(c)). Although this creates a point defect in the reference configuration at which W is infinite, the overall energy is decreased because is much smaller. This is rigorously analyzed in the class of radial deformations in [13]. In fact cavitation is a standard failure mechanism for rubber (see, for example, [77, 114]).

However, maps cannot have planar discontinuities such as occur in cracks. In order to model such fracture we must enlarge the function space again, this time to be the space consisting of special mappings of bounded variation [11]. This space is somewhat technical to describe, but the key points are that (i) can have jump discontinuities at points of a well-defined jump set S_y, at which there is a well-defined unit normal with respect to which there are well-defined limiting values y₊, y₋ from either side of S_y, (ii) the gradient ∇y can still be defined away from S_y. There is, however, an energetic penalty f(y₊, y₋, ν) per unit area associated with a jump in y, so that the energy functional has to be modified to (3.5) for some suitable f, where denotes the element of area on S_y. This functional coincides with I(y) for . Then we expect to obtain a different minimizer with cracks (see Fig. 3(d)). Such ‘free discontinuity’ models for fracture were first proposed by Francfort & Marigo [69]. For a rigorous treatment of fracture in which cavitation is also allowed, together with numerical computations, see [87–90]. In summary, we have an energy functional having different minimizers in three different function spaces.

3.2. The Lavrentiev phenomenon for liquid crystals

But does the Lavrentiev phenomenon occur for liquid crystals? As a first example we can consider the hedgehog defined in (2.11(2.11) ). Since for topological reasons there is no continuous n satisfying on the boundary ∂B of the unit ball B, we have formally that (3.6) whereas .

Lest (3.6(3.6) ) be regarded as somehow cheating, Hardt & Lin [84] constructed an example for the one-constant approximation in which there is a smooth (degree zero) map such that for some α > 0 for all continuous m ∈ H¹(B; S²). In this case there are smooth maps for which the energy is finite, but they all have energy at least α greater than the minimum among H¹(B; S²) maps.

A more interesting possibility is if we allow the director n to be discontinuous across surfaces (planar defects), in the same spirit as the fracture models of Francfort & Marigo, so that we have an augmented functional defined for n ∈ SBV(Ω; S²) by (3.7) for some suitable continuous interfacial energy function f. As before, S_n denotes the jump set of n, its unit normal, and n₊, n₋ the limiting values of n from either side of S_n. Such models have been investigated in [15, 16], where in particular it is shown that invariance requirements imply that f depends only on the quantities and . See also Bedford [28].

To see why the Lavrentiev phenomenon arises for , consider the order reconstruction problem in which a nematic liquid crystal occupies the region Ω = (0, l₁) × (0, l₂) × (0, d) between two parallel plates given by x₃ = 0 and x₃ = d, where d > 0, with antagonistic boundary conditions on the plates (3.8) and periodic boundary conditions on the other faces (3.9) Such problems have been considered by many authors using the Landau - de Gennes model (see Section 4) [23, 24, 32, 40, 110, 144], the variable scalar order parameter model of Ericksen [12], and molecular dynamics [149]. If n ∈ H¹(Ω; S²) satisfies the boundary conditions (3.8(3.8) ), (3.9(3.9) ), then by (2.4(2.4) ) and Jensen’s inequality , where |Ω| = l₁l₂d denotes the volume of Ω, we have that (3.10) However, in the space SBV(Ω; S²) we can consider the competitor (3.11) for which (3.12) which is less than provided (3.13) Thus for sufficiently small plate separation the minimum value for in the space SBV(Ω; S²) is less than that in H¹(Ω; S²). For more details see [15, 16].

The use of SBV(Ω; S²) also has the purely mathematical advantage of treating index defects using a vector field rather than a line-field, by allowing n to jump to − n across surfaces with zero energy cost, such as the surface OP indicated in Fig. 2(c). This is explored in Bedford [28]. In [16] further possible applications to stripe domains in nematic elastomers and smectic thin films are suggested.

4. The Landau - de Gennes model

The Landau - de Gennes model uses a five-dimensional tensor order parameter based on the probability distribution ρ(x, p) of molecular orientations p ∈ S² at a point x. Here p is parallel to the long axis of a molecule, and we regard p and − p as being equivalent. Two generally perceived advantages of the Landau - de Gennes model over that of Oseen-Frank are that (i) it gives structure to defects, so that in particular they have finite energy, (ii) it resolves the problem of orientability of the director.

4.1. The probability distribution of molecular orientations

In order to give meaning to the probability distribution ρ(x, p), consider the ball with centre x and small radius δ > 0. We want δ to be smaller than macroscopic length scales, but large enough to contain many molecules. To give an idea of the numbers, if δ = 1μm then B(x, δ) would contain about a billion molecules. Then ρ(x, p) can be thought of as (a smoothed out version of) the probability of a molecule chosen at random from those in B(x, δ) having orientation p. See Section 5.3 for further discussion.

Since ρ is a probability distribution and ± p are equivalent, ρ satisfies (4.1) where dp denotes the area element on S². By (4.1(4.1) ) the first moment (4.2) The second moment (4.3) is a symmetric second-order tensor which is positive definite, that is M(x)e · e > 0 for all e ∈ S². Indeed (4.4) while M(x)e · e = 0 implies that ρ(x, p) = 0 for p · e ≠ 0, contradicting . Also tr M(x) = 1. The case corresponds to an isotropic molecular distribution at x, and it is easily checked that then .

The de Gennes Q tensor is defined as (4.5) Thus Q(x) measures the deviation of M(x) from its isotropic value, and from the properties of M(x) satisfies (4.6) where λ_min(Q(x)) denotes the minimum eigenvalue of Q(x). Thus Q(x) is a five-dimensional order parameter, whereas the director n(x) is two-dimensional.

4.2. The Landau - de Gennes energy functional

Following de Gennes we suppose that the free energy for a nematic at constant temperature is given by (4.7) Writing H = ∇Q, where Q = (Q_ij), H = (H_ijk) = (Q_{ij, k}), as a consequence of frame-indifference and material symmetry, ψ = ψ(Q, H) should satisfy for any R = (R_ij) ∈ O(3) = {R: R^TR = 1} the isotropy condition (4.8) where Q*_ij = R_irR_jsQ_rs, H*_ijk = R_irR_jsR_ktH_rst.

It is usual to decompose the free-energy density ψ as (4.9) It is often assumed that ψ_B has the form (4.10) where b > 0, c > 0 and a depends linearly on temperature. If then ψ_B is minimized by Q having the uniaxial form (4.11) where (4.12)

Usually it is assumed that ψ_E(Q, ∇Q) is quadratic in ∇Q. Examples of isotropic functions quadratic in ∇Q are the invariants I_i = I_i(Q, ∇Q): (4.13) The first three linearly independent invariants I₁, I₂, I₃ span the possible isotropic quadratic functions of ∇Q. The invariant I₄ is one of 6 possible linearly independent cubic terms that are quadratic in ∇Q (see [30, 126, 150, 158]). Note that (4.14) is a null Lagrangian.

We assume that (4.15) where the L_i are material constants⁴. Note that if L₂ = L₃ = L₄ = 0 then ψ_E = L₁|∇Q|², the one-constant approximation of the Landau - de Gennes theory.

4.3. From Landau - de Gennes to Oseen-Frank

Since ψ_B is minimized for uniaxial Q given by (see (4.11(4.11) )) (4.16) in the limit of small elastic constants L_i (see below for further discussion) we expect minimizers of I(Q) to be nearly uniaxial. This motivates the constrained theory⁵ in which we minimize I(Q) subject to the constraint (4.16(4.16) ) for fixed s > 0. Putting (4.11(4.11) ) into ψ_E we obtain the Oseen-Frank energy with (4.17) The 4 × 4 matrix in (4.17(4.17) ) is invertible, so that the L_i can be also be obtained from the K_i through the formula (4.18) This is one reason why we choose the simplified form (4.15(4.15) ) for the elastic part of the energy, rather than including all the possible cubic terms. Including I₄ as the single cubic term is a common choice (see, for example, [136, 152]), and is convenient for establishing existence when L₄ ≠ 0 (see the discussion at the end of Section 5.2).

As observed by Gartland [74], care needs to be taken as regards the physical interpretation of the L_i being small, because their values depend on the units of length used. A change of units of length can be represented by the change of variables x = λy, λ > 0, , leading to , where (4.19) If the diameter of the physical domain Ω in the x −coordinates is R, and we set λ = R then we can write with . Writing , where e₀ has dimensions of energy per unit volume, (4.19(4.19) ) takes the form (4.20) in which the scaled coefficients are dimensionless. The constrained theory should be applicable when the L′_i are small, which in particular is the case in the large body limit R → ∞ for fixed L_i.

There are a number of rigorous results concerning the passage from the Landau - de Gennes model to that of Oseen-Frank. Most are for the one-constant approximation (4.21) with prescribed uniaxial boundary data, in the limit as L → 0 +, with ψ_B given by (4.10(4.10) ) with a < 0.

Majumdar & Zarnescu [129] showed that for any sequence L^(k) → 0 + there is a subsequence such that minimizers for converge in to a minimizer of subject to the constraint , the convergence meaning that and . If Ω is simply-connected then any uniaxial is orientable [22]. Hence n* is a harmonic map (a minimizer for the one-constant Oseen-Frank energy), and so by Theorem 2 has a finite number of point defects. Away from these defects the convergence is much better. The rate of convergence was improved by Nguyen & Zarnescu [140], who also obtained a first-order correction term.

Canevari [37] studied the convergence of minimizers Q_L as L → 0 + under the logarithmic scaling (4.22) which allows the appearance of line defects in the limit. He shows that these defects consist of straight line segments. (For colloids, where curved disclinations are seen, there are other small geometric parameters. See the study of Saturn rings by Alama, Bronsard & Lamy [8].) In earlier work, Bauman, Park & Phillips [25] studied how line defects appeared in the zero elastic constant limit for thin films, for general elastic constants L₁, L₂, L₃ with L₄ = 0.

Given the appearance of uniaxial Q in the limit L → 0 one can ask whether in fact there are equilibrium solutions for I_L(Q) with L > 0 which are everywhere uniaxial with director n and s = s(x) (the starting assumption of the Ericksen theory of liquid crystals [61]). Lamy [111] shows that in one and two dimensional situations n has to be constant. These results suggest that, although for small L configurations predicted by Landau - de Gennes are very close to being uniaxial (at least away from defects), they are rarely exactly so.

4.4. Description of defects in the Landau - de Gennes model

In contrast to the Oseen-Frank theory, we expect minimizers of I(Q) subject to suitable boundary conditions, such as where is given, to be smooth. When L₄ = 0 this was proved by Davis & Gartland [52] under the conditions (4.23) which imply [126] that ψ_E(Q, ∇Q) ⩾ c|∇Q|² for some c > 0, so that ψ_E(Q, ·) is convex, and the corresponding Euler-Lagrange equations form a semilinear elliptic system.

Thus defects in the Landau - de Gennes theory are not described by singularities in Q. This has been explored for the hedgehog defect by many authors (see, for example, [75, 86, 91, 98, 107, 109, 128]), and for defects of other index in [99, 100].

However L₄ = 0 implies K₁ = K₃, so we need L₄ ≠ 0. But then

Theorem 4

(Ball & Majumdar [18, 19]). If L₄ ≠ 0, and for any continuous ψ_B and any boundary conditions, is unbounded below.

5. Onsager and molecular models

5.1. Onsager models

In the Onsager model, the free energy for a homogeneous nematic liquid crystal at temperature θ > 0 is given in terms of the probability distribution ρ = ρ(p) of molecular orientations by (5.1) As in Section 4, we assume that ρ satisfies the condition of head-to-tail symmetry (5.2) Two well-known examples of suitable kernels k = k(p · q) are

(i)	,
(ii)	,

where κ > 0 is a constant, which may depend on temperature and concentration. For both of these cases the behaviour of solutions depends only on the dimensionless parameter (5.3) The Euler-Lagrange equation corresponding to the free-energy functional (5.1(5.1) ) subject to the constraint is (5.4) where c is a constant Lagrange multiplier, for which one solution valid for all θ > 0 is the isotropic state (5.5) Note that if ρ = ρ(p) is a solution, so is ρ_R(p) = ρ(Rp) for any R ∈ SO(3). In order to study all possible solutions one can begin by looking for solutions bifurcating from . The possible bifurcation points correspond to nonzero solution pairs (ϕ, λ) with to the equation obtained by linearizing (5.4(5.4) ) about , namely (5.6) where λ = −4πk_Bθ. The eigenvalues λ < 0 of (5.6(5.6) ) have been calculated explicitly for a wide range of kernels k by Vollmer [169]. The case of the Maier-Saupe kernel is unusual in that there is only one such eigenvalue , corresponding to in (5.3(5.3) ), where there is a transcritical bifurcation, and for this kernel all critical points can be explicitly determined (Fatkullin & Slastikov [66], Liu, Zhang & Zhang [125]).

For the Onsager kernel, there are infinitely many distinct bifurcation points α_j. Vollmer [169] (see also Wachsmuth [170]) shows that at the least bifurcation point there is a transcritical bifurcation to an axially symmetric state, together with rotations of it, and proves other properties of the solution set. However a complete description of the set of all solutions remains an open problem. For the 2D case see Lucia & Vukadinovic [127].

5.2. The singular bulk potential

The Q tensor corresponding to ρ is (5.7) Following Katriel, Kventsel, Luckhurst & Sluckin [103] and Ball & Majumdar [18, 19], the bulk energy ψ_B(Q, θ) can be identified with the minimum of I_θ(ρ) among probability distributions ρ on S² satisfying (5.2(5.2) ) such that Q(ρ) = Q.

For the Maier-Saupe kernel the interaction term takes the simple form (5.8) so that (5.9) Hence (5.10) where (5.11)

Theorem 5

([18, 19]). f is a smooth strictly convex function of Q defined for , and (5.12) for constants C₁, C₂.

In particular the singular potential ψ^s_B(Q, θ) = ψ_B(Q, θ) satisfies (5.13) Using this bulk potential, under suitable conditions on the L_i it follows that the energy I(Q) in (4.7(4.7) ) is bounded below and attains a minimum even for L₄ ≠ 0, in contrast to Theorem 4. To see this, note that implies that (5.14) Define (5.15) Thus by (4.23(4.23) ) if (5.16) we have that ψ_E(Q, ∇Q) ⩾ c₁|∇Q|² for some c₁ > 0, implying in particular that the K_i given in (4.17(4.17) ) satisfy the Ericksen inequalities (2.3(2.3) ), as can be checked directly.

Furthermore the singular potential predicts an isotropic to nematic phase transition in the same way as the quartic potential (4.10(4.10) ), which can be regarded as an expansion of the singular potential around the isotropic state Q = 0.

For a general framework for singular potentials together with applications see the recent work of Taylor [167].

5.3. Molecular models

It is an important open problem to establish some kind of rigorous link between molecular models of liquid crystals, however primitive, and continuum models, so that solutions to the equations of motion for the molecules are shown to converge in an appropriate average sense to those for some continuum model. This is a deep problem even for the much simpler situation of a moderately rarified gas, for which there is a large literature (see, for example, [70, 112]) on the passage from Newtonian mechanics of hard spheres to a description of the dynamics of the gas in terms of a probability distribution for the velocities satisfying the Boltzmann equation, and a further large literature (see, for example, [78, 155]) on how to pass from the Boltzmann equation to the equations of continuum fluid dynamics. There do not seem to be any corresponding rigorous results known for liquid crystals, though there are penetrating formal studies (for example, [76, 139]).

One issue is how the macroscopic variables, such as ρ(p) and its moments, would emerge in such a theoretical framework. We already gave in Section 4.1 a suggestive description of how this might happen for nematics via a coarse-graining procedure. Less clear is the situation for smectics, in some theories of which (for example [44, 54, 79, 105, 133, 148]) the molecular mass density c(x) plays a role, its variation describing smectic layers. How can we have a macroscopic variable which varies on a microscopic length scale? One might try to resolve this by averaging at a fixed time over a ball centred at x and small radius δ sufficient to contain lots of molecules. However, it is not obvious this works because such an average may not detect oscillations. For example, (5.17) is zero for all x if kδ is a root of the equation sin λ = λcos λ. This is related to the Pompeiu problem [7, 151], which asks for which bounded domains is it true that the only locally integrable function for which is f = 0, for which (5.17(5.17) ) is a counterexample for Ω a ball. A possible remedy would be to average also over a short time interval, so that the radius δ could be made smaller than the layer thickness. To make this rigorous would seem to require a quite detailed understanding of the molecular motion in smectics.

6. Unequal elastic constants

We have seen that much rigorous mathematical work on liquid crystals makes the one-constant approximation (both for Oseen-Frank and Landau - de Gennes). In this section we highlight some rigorous results for the case of unequal elastic constants.

6.1. Oseen-Frank

One reason why the case of general K_i may be more difficult is that it is only for the one-constant approximation that the Lagrange multiplier λ(x) in the Euler-Lagrange equation (2.8(2.8) ), that is depends only on ∇n.

In the one-constant approximation (Theorem 2) energy minimizers are smooth except for a finite number of point defects. This is not known for general K_i. However, we have

Theorem 6

(Hardt, Lin & Kinderlehrer [81]). Let K₁ > 0, K₂ > 0, K₃ > 0, and let . Then any minimizer n ∈ H¹(Ω; S²) of subject to , where W is given by (2.2(2.2) ), is analytic outside a closed subset of Ω whose Hausdorff dimension is less than one.

Hardt, Lin & Kinderlehrer [81] also prove a corresponding partial regularity result up to the boundary.

An extension of Rivière’s result for the one-constant approximation to general K_i is due to Hong [93], who proves that if and if ∂Ω is smooth with smooth on ∂Ω and nonconstant, then there are infinitely many solutions to (WEL) satisfying (see Section 6.1.1 for the special case of the hedgehog on a ball).

6.1.1. Energy minimizing properties of universal solutions

It is interesting to ask under which conditions on the K_i, and for which boundary conditions, the universal solutions described in Section 2.4 are energy minimizers. For the hedgehog , it was observed by Hélein [85] that the proof of Lin [118] works in the case K₂ ⩾ K₁ > 0, K₃ > 0, so that under these conditions on the K_i the hedgehog minimizes I(n) subject to its own boundary conditions. A detailed proof of this for Ω a ball is given in [141]; for a general bounded domain Ω one can extend any n ∈ H¹(Ω; S²) with to a larger ball B by setting outside Ω, and then apply the result for B, but a direct proof is also possible. The hedgehog is a pure splay configuration, since . The condition that K₂ ⩾ K₁ says that it is energetically easier to splay than to twist.

Necessary and sufficient conditions for to be a minimizer subject to its own boundary conditions do not seem to be known. Hélein [85] showed that if 8(K₂ − K₁) + K₃ < 0 then the second (outer) variation of I(n) at can be negative, so that is not a minimizer. On the other hand Cohen & Taylor [49] showed that if the opposite inequality 8(K₂ − K₁) + K₃ ⩾ 0 holds then the second variation of I(n) at is strictly positive, so that with respect to certain variations is a local minimizer. A simplified proof of these second variation calculations was given by Kinderlehrer & Ou [104]. Alouges & Ghidaglia [10] presented numerical computations indicating that the condition 8(K₂ − K₁) + K₃ ⩾ 0 is not sufficient for to be a global minimizer; however, K₁ = K₂ for two of the three cases when a configuration with lower energy was apparently found, when by the result mentioned above is in fact a minimizer.

Hong [92] shows that for a ball B there are always at least two solutions of (WEL) satisfying that are partially regular (that is, smooth outside a closed subset of B of Hausdorff dimension less than 1), while if K₁ ⩽ min (K₂, K₃) there are infinitely many such partially regular solutions having non-negative second variation.

One can similarly study when pure twists of the form (2.10(2.10) ) are minimizers. As in Section 3.2 consider the region Ω = (0, l₁) × (0, l₂) × (0, d) between two parallel plates given by x₃ = 0 and x₃ = d > 0, with planar boundary conditions on the plates (6.1) where n₀ = (cos α, sin α, 0), n_d = (cos β, sin β, 0) and α, β ∈ [0, π), and periodic boundary conditions (6.2) on the other faces, where we have expressed the boundary conditions in terms of n⊗n because ± n are physically indistinguishable (for simplicity we didn’t do this in Section 3.2). Then we have

Theorem 7.

The minimum of I(n) among n ∈ H¹(Ω; S²) satisfying the boundary conditions (6.1(6.1) ), (6.2(6.2) ) is attained. If K₂ ⩽ min (K₁, K₃) then any minimizer is a pure twist of the form (6.3) where (6.4)

Thus up to a physically meaningless change of sign, there is a unique minimizer if n₀ · n_d ≠ 0, and two minimizers if n₀ · n_d = 0 (corresponding to twists in opposite directions). The proof, which will appear elsewhere, is a consequence of the identity (2.5(2.5) ) and the fact that for n satisfying the boundary conditions The conclusion of the theorem is intuitively reasonable, since the condition K₂ ⩽ min (K₁, K₃) says that it is energetically easier to twist than to splay or bend.

If K₂ > min (K₁, K₃) then pure twists are not in general minimizers. For example, if K₁ < K₂ = K₃ a twist-bend equilibrium solution depending only on x₃ has less energy than any pure twist (see Leslie [117], Stewart [165]). Using a deflation (see [65]) numerical scheme Adler et al. [5] compute a variety of equilibrium solutions for some different cases, studying for example the bifurcation to twist-bend solutions as K₂ is increased for fixed K₁, K₃. It would be interesting to understand why (the presumably infinitely many) less regular equilibria similar to those whose existence is proved by Rivière and Hong are not found by the implementation of the deflation method.

6.2. Landau - de Gennes

Consider the problem of minimizing the free energy (4.7(4.7) ) with the singular bulk potential subject to , where is smooth and for all x ∈ ∂Ω.

As explained at the end of Section 5.2, if the inequalities (5.16(5.16) ) hold, then the minimum of I(Q) is attained. Let Q* be a minimizer, so that in particular we have that (6.5) Is Q* smooth? This would be straightforward to prove if we knew that (6.6) for some constant ϵ > 0. For then we could show that Q* is a weak solution of the corresponding Euler-Lagrange equations and use elliptic regularity theory in the same way as Davis & Gartland [52].

But surely (6.6(6.6) ) must be true, because why would it be good for the integrand in (6.5(6.5) ) to be infinite somewhere, when Q* a minimizer? However in fact this phenomenon often arises in the calculus of variations, as we have already seen for the hedgehog defect, for the model of cavitation discussed in Section 3, and more surprisingly in the one-dimensional example mentioned there, so it is a delicate matter to show that it cannot happen.

A proof of (6.6(6.6) ) in the one-constant approximation ψ_E(Q, ∇Q) = L₁|∇Q|² is given in [18, 19], but the method seems not to work for more general elastic constants, for which the question of whether (6.6(6.6) ) holds remains open. However, without proving (6.6(6.6) ), Evans, Kneuss & Tran [64] show as a consequence of a more general partial regularity result that Q* is smooth outside a closed subset E⊂Ω of measure zero. Bauman & Phillips [26] study the regularity problem in 2D, in particular proving (6.6(6.6) ) when L₂ = L₃ = 0 and L′₁ > 0, where L′₁ is defined in (5.15(5.15) ).

7. Omissions

This paper does not pretend to be a comprehensive survey of mathematical work on liquid crystals. Some notable omissions concern:

1.	Dynamics: Here one goal is to understand the qualitative properties of solutions to the Ericksen-Leslie equations, which were briefly mentioned in Section 2.6. These equations consist of a momentum equation generalizing the Navier-Stokes equations for motion of a linear viscous fluid, coupled to an evolution equation for the director n. This system of equations poses additional challenges to those already present for the Navier-Stokes equations, for which it is famously not known whether solutions to initial-boundary value problems in 3D are smooth or unique. A particular added difficulty for the Ericksen-Leslie equations is how to handle the unit vector constraint |n| = 1. Thus the first attempts to prove existence of solutions relaxed this constraint by replacing W(n, ∇n) by W(n, ∇n) + ϵ− 2(|n|2 − 1)2 in the hope of recovering the constraint in the limit ϵ → 0 +, and global existence and partial regularity of weak solutions to a simplified version of the corresponding system of equations in 3D was proved in [119, 121, 122]. For versions of the Ericksen-Leslie equations in 2D with the unit vector constraint there is a relatively complete global existence, uniqueness and regularity theory (see, for example, [94–96, 120, 123]). These papers ignore the moment of inertia of molecules; a short-time existence theory when this is included but without dissipation in the director equation is provided in [43]. For special initial data global existence of a weak solution in 3D to a simplified Ericksen-Leslie system was proved in [124], and for the same system the existence of solutions developing a singularity in finite time was proved in [97]. These papers mostly concern the one-constant approximation. For a discussion of numerical approximation of solutions to the Ericksen-Leslie equations see [171]. There are corresponding sets of dynamical equations when the order paranmeter is the Q tensor, such as the Beris-Edwards model [29]. For recent work on the existence and properties of solutions see for example [1, 2, 142, 143], and in the context of the singular bulk potential [67, 174]. For results on how to relate these models to the Ericksen-Leslie equations see [172]. The dynamical behaviour of liquid crystals has many complex aspects (see e.g. [68]) worthy of mathematical treatment.
2.	Smectics. The results highlighted in this paper mainly concern nematics. Among mathematical studies of other liquid crystal phases, smectics (briefly touched on in Sections 3.2, 5.3) have generated considerable recent interest. See, for example, [9, 27, 31, 36, 45–47, 50, 71–73, 79, 101, 102, 134, 145–147, 156, 162].
3.	Liquid crystal elastomers. Liquid crystal elastomers are materials formed from polymers to which are attached liquid crystal mesogens. In the case of nematic elastomers the model of Bladon, Warner & Terentjev [33, 173] is based on a free-energy density ψ(∇y, n) depending on the deformation gradient ∇y and n, but not on ∇n (so that Frank elasticity is ignored). Minimizing with respect to n leads to a free-energy density W(∇y) = min |n| = 1ψ(∇y, n) which was shown by De Simone & Dolzmann [56] not to be quasiconvex, leading to nonattainment of energy minimizers and an explanation of observed laminated structures [108] similar to those seen in martensitic phase transformations. For a selection of recent mathematical work in this area see [6, 35, 41, 42].
4.	Topological aspects. Liquid crystals traditionally have been an area in which geometry and mechanics meet topology, for example in the topological description of defects. An interesting area of current activity in which topology plays a role is that of nematic shells (see, for example, [38, 39, 154, 160, 161]). Related to this work is the interest in the topology of disclination lines induced by colloidal systems (see, for example, [8, 51, 132, 166, 168]). Other interesting connections with topology appear in [3, 4].

Acknowledgments

I acknowledge with thanks the hospitality of the Liquid Crystal Research Institute at Kent State University, where much of this lecture was prepared. I am grateful to Giacomo Canevari, Patrick Farrell, David Kinderlehrer, Fanghua Lin, Tom Lubensky, Peter Palffy-Muhoray, Dan Phillips, Michaela Vollmer and Arghir Zarnescu for illuminating discussions and references. I especially thank Epifanio Virga and the referee for their careful reading of the article and very useful suggestions and comments.

Funding

This research was supported by EPSRC (GRlJ03466, the Science and Innovation award to the Oxford Centre for Nonlinear PDE EP/E035027/1, and EP/J014494/1), the European Research Council under the European Union’s Seventh Framework Programme (FP7/2007-2013) / ERC grant agreement no 291053 and by a Royal Society Wolfson Research Merit Award.

Notes

1 In MathSciNet, the main mathematics reviewing journal, the number of papers mentioning liquid crystals rose from 4 for 1971–1975, to 223 for 2006–2010, and 407 for 2011–2015.

2 I am grateful for discussions with Fanghua Lin and Epifanio Virga conerning the material in this section.

3 for example, the problem of minimizing the integral I(u) = ∫¹_{− 1}[(x⁴ − u⁶)²u_x²⁸ + ϵu²_x] dx subject to the boundary conditions u( − 1) = −1, u(1) = 1, where ϵ > 0 is sufficiently small. For a detailed discussion, together with other historical references, see [20].

4 Here we use the definitions of the L_i, i = 1, 2, 3 most common in the literature (for example [53, 126, 135]), rather than those used in the papers [16, 19, 21, 22] in which the L_i were permuted with respect to these definitions, with a corresponding permutation of the I_i.

5 For a corresponding constrained theory for biaxial nematics see [137, 138].