I am trying to understand what the mean field approximation means when expressed in tensor notation for the Maier-Saupe Model of nematic liquid crystals. I am following along with Jonathan Selinger's lecture and my question appears around 16 mins in. The system in question is a system of cylindric molecules orientationaly aligned along a director. My confusion arises from how the angular dependence of this potential is expressed and how the mean field approximation is made concerning this angle.

To begin with, we have pairwise interaction energy: $$ V_{int} = -J P_2(\cos(\gamma)) = -J(\frac{3}{2}\cos^2\gamma - \frac{1}{2})$$

where J is the coupling strength of the interaction and $\gamma$ is the angle between the cylindrical molecules.

This pairwise interaction energy is an anisotropic potential which pushes the molecules to align with one another. Given this interaction potential, we want to calculate the free energy of the system with respect to an order parameter to see if the system is isotropic or nematic for a given temperature.

To calculate the free energy $F = \langle E \rangle - TS$, we take the average of the internal energy

$$ \langle E \rangle = - J \frac{Nq}{2} \langle \frac{3}{2}\cos^2\gamma - \frac{1}{2} \rangle$$

Where $\langle E \rangle $ denotes the average internal energy and $\frac{Nq}{2}$ is the number of interacting pairs. The cosine squared term is then expressed as a dot product of the two unit vectors of each molecule, $\hat{l}$ and $\hat{m}$, and that dot product is then expressed in Einstein summation notation for tensors, which is still fine with me:

$$ \langle \frac{3}{2}\cos^2\gamma - \frac{1}{2} \rangle = \langle \frac{3}{2}(\hat{l} \cdot \hat{m})^2 - \frac{1}{2} \rangle = [ \frac{3}{2} \langle l_im_il_jm_j \rangle - \frac{1}{2}]$$.

My issue is with the next step. Now, Selinger breaks up this average of products into a product of averages. However, it is broken so that only terms from one molecule are multiplied by another:

$$\langle l_im_il_jm_j \rangle = \langle l_il_j \rangle \langle m_im_j\rangle $$

So he says here that this is the mean field approximation. That mathematically this is equivalent to assuming the molecules are uncorrelated with one another. Then once breaks things up into these outer products, it is easy to relate these outer products to the nematic order parameter of the Q tensor, and do the desired analysis.

My questions are:

  1. How does this final step correspond to making a mean field approximation?
  2. If this is just an approximate formula, what would the full formula be for this step? I would really appreciate being able to to explicitly see what terms we are neglecting with this approximation.
  3. If I understand correctly,$\langle l_il_j \rangle \langle m_im_j \rangle$ just implies regular old matrix multiplication here to generate a 3x3 tensor. How are we taking a scalar (the cosine squared of the intermolecular angle) and approximating this with a 3x3 tensor?

The explanation Selinger gives for this step is related to the mean field approximation of the Ising model, where he draws a connection between the approximation in the Ising model, where $\pm 1 = \sigma_i$ for the spin at site i, and the mean field approximation is expressed as $\langle \sigma_i \sigma_j \rangle = \langle \sigma_i \rangle \langle \sigma_j \rangle$.

These notes draw a connection between this approximation and the covariance/Pearson correlation coefficient where $\rho(X,Y) = \frac{Cov(X,Y)}{\sigma_x \sigma_y} = \frac{\langle XY \rangle - \langle X \rangle \langle Y \rangle}{\sigma_x \sigma_y} $ where $\sigma_x$ is now the std dev of x and $\rho(X,Y)$ the pearson correlation coefficient. He points out if $\rho(X,Y) = 0$ the two variables are said to be "uncorrelated" which implies $ \langle XY \rangle - \langle X \rangle \langle Y \rangle = 0 \rightarrow \langle XY \rangle = \langle X \rangle \langle Y \rangle$, which has the form of the approximation we are making. However, even given this level of explanation, it is still not clear to me (1) how this formula relates to the tensor manipulations above and (2) how saying the two molecule's orientations are uncorrelated corresponds to the physical intuition behind the mean field approximation, that being that you are "averaging out" all the interactions felt by one molecule into a single field of interaction for that molecule. Wouldn't this average field still cause some correlation in the orientations of the molecules? Indeed isn't that the point of this method - to generate an average potential which causes the orientations to be correlated with one another along a director to generate a nematic phase?

  • 1
    $\begingroup$ Use \langle and \rangle for $\langle$ and $\rangle$, respectively. $\endgroup$ Sep 27, 2023 at 8:56

1 Answer 1


My attempt to answer the questions.

As I understand it, the classical system of molecules is being considered here. The state of each molecule is described by a three-dimensional vector of unit length. The state of a system of $N$ molecules is described by $N$ vectors $\vec{l^1},\ldots,\vec{l^N}$. The aim is to investigate the emergence of an orientation order in such a system. Therefore, the following model simplifications are used. The molecules are considered to be placed in the nodes of a regular lattice, and only the interaction of the nearest neighbors is taken into account. According to the chosen pair interaction, the energy of a particular configuration of a system of molecules is $$ E_N(\vec{l^1},\ldots,\vec{l^N}) = -J\sum_{<a,b>}\left(\frac32(\vec{l^a},\vec{l^b})^2 - \frac12\right) \quad(1) $$ Vectors $\vec{l}$ are continuous variables, so we have to consider the probability density of configurations of molecules. If the system is at a temperature $T$, then this probability density has the form of a Gibbs distribution $$ W_N(\vec{l^1},\ldots,\vec{l^N}) = C\exp\left(-\frac1{k_bT} E_N(\vec{l^1},\ldots,\vec{l^N}) \right) = $$ $$ = C\ \exp\left(\frac{J}{k_BT} \sum_{<a,b>}\ \left(\frac32(\vec{l^a},\vec{l^b})^2 - \frac12\right) \right),\quad (2) $$ where $C$ is the normalization constant. In statistical mechanics, the internal energy of a system is defined as the average value of the configuration energy (1) according to the Gibbs distribution (2) $$ \langle E \rangle = \int\!\ldots\!\int\ \prod_{с = 1}^{N} d\vec{l^с}\ E_N(\vec{l^1},\ldots,\vec{l^N}) W_N(\vec{l^1},\ldots,\vec{l^N}) = $$ $$ = -J\sum_{<a,b>}\ \int\!\ldots\!\int\ \prod_{с = 1}^{N} d\vec{l^с}\ \left(\frac32(\vec{l^a},\vec{l^b})^2 - \frac12\right)\ W_N(\vec{l^1},\ldots,\vec{l^N}), \quad (3) $$ It is impossible to calculate the internal energy (3) and even the normalizing constant $C$ in the case of distribution (2). Therefore, we need to use approximations, in particular, the mean field theory.

In the mean field theory, the Gibbs distribution for a system of a large number of interacting particles is approximated by the distribution for a system of independent particles. That is, when calculating the internal energy (3), the exact distribution (2) is replaced by an approximate: $$ W_N^{mf}(\vec{l^1},\ldots,\vec{l^N}) = \prod_{a=1}^N \rho(\vec{l^a}),\quad (4) $$ where the one-particle distribution $\rho(\vec{l})$ is defined in a self-consistent manner using the mean field equations. According to distribution (4), $\vec{l^a}$ and $\vec{l^b}$ for $a \neq b$ are identically distributed independent random variables. In particular, this means that the average value of any product of the functions of the vectors $\vec{l^a}$ and $\vec{l^b}$, $a\neq b$, is equal to the product of the average values of these functions: $$ \langle f(\vec{l^a}) g(\vec{l^b}) \rangle = \int\!\ldots\!\int\ \prod_{с = 1}^{N} d\vec{l^с}\ f(\vec{l^a}) g(\vec{l^b}) W_N^{mf}(\vec{l^1},\ldots,\vec{l^N}) = $$ $$ = \int d\vec{l^a}\ f(\vec{l^a})\rho(\vec{l^a})\ \int d\vec{l^b}\ g(\vec{l^b})\rho(\vec{l^b}) = \langle f(\vec{l^a}) \rangle \langle g(\vec{l^b}) \rangle %\quad (5) $$ Therefore the averaging of the quantity $(\vec{l^a},\vec{l^b})^2 = \sum_{i,j=1}^3 l^a_i l^b_i l^a_j l^b_j$ over the distribution (4) gives $$ \langle (\vec{l^a},\vec{l^b})^2 \rangle = \sum_{i,j=1}^3 \langle l^a_i l^a_j \rangle \langle l^b_i l^b_j \rangle = \sum_{i,j=1}^3 Q_{ij} Q_{ij} = \mbox{Tr}\ Q^2, $$ where $Q$ is a $3\times 3$ matrix with elements $$ \langle l^a_i l^a_j \rangle = \langle l^b_i l^b_j \rangle = Q_{ij} = \int\! d\vec{l}\ l_il_j\ \rho(\vec{l}). $$ So we not only have regular old matrix multiplication here, but also the operation of taking a trace. The approximate expression for the internal energy is now $$ \langle E \rangle \approx -J\frac{Nq}2 \left(\frac32 \mbox{Tr}\ Q^2 - \frac12 \right) $$

In the mean field theory, the optimal method for determining the approximate distribution (3) is known, associated with minimizing the upper bound of Bogolyubov (Feynman) for free energy. Nevertheless, it is difficult to quantify exactly what we neglected when we replaced the exact distribution (2) with the approximate distribution (3). To do this, we actually need to find an exact expression for the internal energy $\langle E \rangle$, which is impossible. Possible way to find approximate estimate of the mistake of mean field method is to write few terms of perturbation expansion on $\log(W_N/W_N^{mf})$.

  • $\begingroup$ What function do the integrals have in (1)? Also, is it correct to say that $W_N()$ is the probability of finding the system in the state specified by the integrals? Additionally, how does it follow that given the fact that the two unit vectors are independently distributed random variables, their averaging can be broken up in that way. Thank you for your answer as a whole, especially the clear explanation to my 3rd original question about the trace! $\endgroup$
    – McKinley
    Sep 28, 2023 at 22:48
  • 1
    $\begingroup$ @McKinley I edited my answer. $\endgroup$
    – Gec
    Sep 30, 2023 at 9:51
  • $\begingroup$ I follow your explanation that $<f(\vec{l^a})g(\vec{l^b})> = <f(\vec{l^a})><g(\vec{l^b})> $. My remaining issue is that $<(\vec{l^a} \cdot \vec{l^b})^2>$ is not an average of a product of two functions (like $<f(\vec{l^a})g(\vec{l^b})>$), but rather an average of a function of two variables: $<f(\vec{l^a} , \vec{l^b})>$. And it seems like derivation is saying $<f(\vec{l^a} , \vec{l^b})> = <f(\vec{l^a})><g(\vec{l^b})>$. How do you reconcile this? Am I wrong in seeing the dot product as a single function of two variables? Thank you so much, I feel like I am so close to understanding. @Gec $\endgroup$
    – McKinley
    Oct 2, 2023 at 22:40
  • 1
    $\begingroup$ @McKinley Square of dot product, $(\vec{l^a},\vec{l^b})^2$, cannot be written as a product of two functions, but it is a sum of product terms, $\sum_{i,j} l^a_i l^a_j l^b_il^b_j$. For each term you have $\langle l^a_i l^a_j l^b_il^b_j \rangle = \langle l^a_i l^a_j\rangle \langle l^b_il^b_j\rangle$. And the average of the sum is the sum of averages, $\langle (\vec{l^a},\vec{l^b})^2 \rangle = \sum_{i,j} \langle l^a_i l^a_j\rangle \langle l^b_il^b_j\rangle$. $\endgroup$
    – Gec
    Oct 3, 2023 at 5:51
  • 1
    $\begingroup$ To clarify for others, $<\sum_{i,j} l^a_il^b_il^a_jl^b_j> = <\sum_{i,j} f(\vec{l^a})g(\vec{l^b})> = \sum_{i,j} <f(\vec{l^a})><g(\vec{l^b})>$. @Gec thank you for all the time and effort you put into answering my questions! $\endgroup$
    – McKinley
    Oct 3, 2023 at 15:57

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.