Maier-Saupe nematogenic fluid with isotropic Yukawa repulsion at a hard wall: Mean field approximation

The mean field approximation is formulated within the framework of the density field theory to study the properties of a Maier-Saupe nematogenic fluid near a hard wall. The density and the order parameter profiles are obtained using the analytical expressions derived in the linearized mean field approximation. The temperature dependencies of the contact values of the density and order parameter profiles are analyzed in detail. To estimate a validity of the applied approximations, the obtained theoretical results are compared with the original computer simulation data.


Introduction
It is a great pleasure and a big honor for us to contribute this paper to the festschrift dedicated to Professor Stefan Sokołowski. Stefan is a well-known expert on the modelling of physico-chemical properties of complex fluids such as chemically reacting fluids in disordered porous media [1,2], studies of the isotropic-nematic phase transition in confined nematic fluids [3], modelling of the properties of fluids in pores with walls decorated with tethered polymer brushes [4], treatment of phase behavior in confined functional colloids [5] and many other complex fluids in complex confinements.
In this paper, we study the influence of the surface on a nematic fluid near a hard wall. Due to orientational ordering, such systems have many unique properties which are very important in the display industry. The anchoring phenomena are among them, according to which the surface induces a specific orientation of the nematic director with respect to the surface [6]. In order to understand the connection between the anchoring phenomena and the interaction between nematic molecules and the surface, there the Henderson-Abraham-Barker (HAB) approach previously developed for isotropic fluids at the wall [7] has been employed [8,9]. In this approach, the distribution of the fluid near the wall is described by the wall-particle Ornstein-Zernike (OZ) equation with the fluid distribution function in the bulk calculated in the mean spherical approximation (MSA) [10]. Using this approach, it is possible to evaluate the role of an orientation-dependent interaction of nematic molecules with the surface in the anchoring phenomena. However, in the MSA, the HAB approach does not correctly take into account the contribution from long-range intermolecular interactions. As a result, it does not satisfy an exact relation known as the contact theorem, which was formulated in [11,12] for isotropic fluids near a wall, and recently in [13] it was reformulated for anisotropic fluids near a hard wall. According to this theorem, the contact value of the density profile of a fluid near a hard wall is determined by the pressure of the fluid in the bulk. In [13], the contact theorem was also formulated for the order parameter profile.
An alternative way of describing fluids at a hard wall was developed within the framework of the density field theory. In this theory, the contributions from the mean field and from fluctuations are separated. The theory was successfully applied to ionic fluids at a hard wall [14][15][16][17] and to simple fluids with Yukawa-type interactions near a hard wall [18,19]. It was shown that the mean field treatment of a Yukawa fluid at a hard wall reduces to the solution of a non-linear differential equation for the density profile, while the treatment of fluctuations reduces to the OZ equation with the Riemann boundary condition [20]. The density field theory was applied to the description of bulk properties of a nematic fluid in [21,22]. The application of the density field theory to the description of a nematic fluid at a hard wall was initiated in [23].
In this paper, we use the mean field approximation to investigate the effect of a hard wall on the properties of a nematic fluid. We demonstrate the principal difference between the behavior of the order parameter profile obtained in the mean field approximation and its linearized version. In order to check the validity of both approaches we compare the obtained theoretical results with computer simulations data.
It is necessary, in numerical calculations, to cut-off the potential ν(r 12 , Ω 1 Ω 2 ) at some finite distance and due to this in expression (1), ν 0 (r ) and ν 2 (r ) are replaced byν i (r ) = ν i (r ) for r r c andν i (r ) = 0 for r > r c , where i = 0, 2 and r c is the cut-off radius.
Within the field-theoretical formalism, the Hamiltonian is a functional of the density field and can be written as a sum of the entropic and the interaction terms 12 , Ω 1 Ω 2 )ρ(r 1 , Ω 1 )ρ(r 2 , Ω 2 )dr 1 dr 2 dΩ 1 dΩ 2 , where β = 1/k B T is the inverse temperature, dΩ = (1/4π) sin θdθdφ is the normalized angle element, ρ(r, Ω) is particle density per angle such that ρ(r, Ω)dΩ = ρ(r), Λ T is the thermal de Broglie wavelength of the molecules, the quantity Λ −1 R is the rotational partition function for a single molecule [26].

Mean field approximation
In this paper, we restrict our consideration to the mean field (MF) approximation which is the lowest order approximation for the partition function. In the canonical formalism, it corresponds to fixing the Lagrange parameter λ such that the following relation is true for the singlet distribution function As a result, we have ρ(r 1 , Ω 1 ) = ρ bulk (Ω 1 ) exp −β ν(r 12 , Ω 1 Ω 2 ) ρ(r 2 , Ω 2 ) − ρ bulk (Ω 2 ) dr 2 dΩ 2 , where is the singlet distribution function for the bulk nematic fluid in the MF approximation, defined within the framework of the Maier-Saupe theory [24,25], κ 2 2 = 4πρ b βA 2 , ρ b is the bulk value of the fluid density, S b = (1/ρ b ) 1 0 P 2 (cos θ)ρ bulk (Ω)d cos θ is the bulk value of the orientational order parameter. After integration with respect to orientation Ω 2 , we obtain ρ(r 1 , Ω 1n , Ω w n ) where Ω w n denotes the angle between the nematic director and the surface, and the mean field potentials V 0 (r 1 , Ω w n ) = β ν 0 (r 12 )ρ(r 2 , Ω w n )dr 2 , V 2m (r 1 , Ω w n ) = β ν 2 (r 12 )S 2m (r 2 , Ω w n )dr 2 .
The bulk values of these potentials are V b is the density profile. The quantities where S * 2m (r, Ω w n ) are the order parameter profiles. Far from the wall we have S * 20 (r, Ω w n ) → S b , S * 2m (r, Ω w n ) → 0 for m 0. Simple calculations show that, in order to take into account the cut-off radius r c , one should substitute the quantities κ 2 i byκ 2 i such that

Linearized MF approximation
The gradient of equation (6) gives where we define an equivalent of the electric field as According to the properties of the Yukawa potential we can write △ − α 2 2 V 2m (r, Ω w n ) = −4πβA 2 S 2m (r, Ω w n ). (15) Due to translational invariance parallel to the wall, the functions considered depend only on the distance z to the wall. Equations (10)-(15) make a set of six differential equations for the unknown functions 13608-3 ρ(r, Ω 1n , Ω w n ), S 2m (r, Ω w n ), E 0 (r, Ω w n ), E 2m (r, Ω w n ), V 0 (r, Ω w n ), V 2m (r, Ω w n ). We note that in the case when the director is oriented perpendicularly to the wall, Ω w n = 0, the singlet distribution function is axially symmetric. Consequently, the equations considered will retain only the terms with m = 0. In this paper, we will restrict our further investigation to this special case.
The resulting solutions of the linearized profile are as follows: where Parameters λ 0 and λ 2 are identical to the parameters found in the bulk phase when Gaussian fluctuations are taken into account [22] and characterize a decay of the isotropic repulsive and the anisotropic attractive interactions, respectively.
Hereafter, the approach presented in this subsection is referred to as the linearized mean field (LMF) approximation. The expressions for ρ(z) and S * 20 (z) obtained within this approximation correspond to the case of infinite cut-off radius r c → ∞.

Contact theorem
As it was shown in [13], the density and order parameter profiles satisfy some exact relations known as the contact theorems. According to these relations in the absence of wall-particle interactions, the contact values of the density profile ρ(z = 0) and of the order parameter profile S 20 (z = 0) do not depend on the angle Ω w n and are equal to ρ(z = 0) = β dΩ 1n P (Ω 1n ) = βP, (22) S 20 (z = 0) = β dΩ 1n P 2 (cos θ)P (Ω 1n ), (23) where P is the bulk pressure and P (Ω 1n ) can be treated as the bulk partial pressure for molecules with a given orientation Ω 1n .
In the MF approximation for the model under consideration, relations (22) and (23) give These relations are used as boundary conditions in the solution of differential equations of the LMF approximation. In order to take into account the cutoff distance r c in relations (24)-(25), we should change κ 2 i toκ 2 i given by equation (11). However, we should note the principal difference between the exact results (22)- (23) and results (24)-(25) of the MF approximation. In paper [21], an invariant was found which was used to prove the contact theorem for the MF density profile of a nematogenic fluid in form (24). However, no similar proof of the contact theorem for the MF order parameter profile in form (25) exists.

Numerical calculation details
In order to verify the theoretical approaches presented in the previous section, a series of numerical calculations were carried out. To this end, for the model pair potential (1), the following parameters were chosen: A 0 /|A 2 | = 3.0 and α 0 /α 2 = 1.6. It should be noted that all quantities marked by a star in our paper are considered as non-dimensional. For instance, all distances are reduced as r * = r α 2 or z * = zα 2 , densities are reduced as ρ * = ρ/α 3 2 and temperature as T * = k B T /(|A 2 |α 2 ). The potential (1) is characterized by a rather soft repulsive part (isotropic contribution) and a small attractive part dependent on the relative orientation between a pair of fluid particles (anisotropic contribution). It is worth noting that the considered fluid is mostly a repulsive one, and a small attractive contribution affects mainly an orientational properties of the fluid, at least at the conditions used in our study. In particular, we consider the fluid at the density ρ * b = 1.0, and the temperature interval T * = 0.5−3.5 is chosen. We have found that, at these temperatures, the considered fluid is beyond its vapour-liquid phase transition region. Therefore, only a nematic-isotropic phase transition can be expected in our case.
To obtain a numerical solution of the integral equation (6) used in the MF approximation, the Picard iterative method was applied. The problem was considered in the cylindrical coordinates, which are set along z-axis normal to the wall surface. The integrations over z were performed using trapezoidal rule with a step ∆z = 0.02/α 2 , while all integrations over r were done analytically. A step of integration over cos(θ) was chosen as 0.0025. To take into account the confinement, we consider a fluid between two hard walls at a distance L z = 36/α 2 to each other. Two cut-off radii r c = 6.0/α 2 and 12.0/α 2 are used in our study. The chosen distance L z appears to be sufficient to get a bulk-like region of a fluid in the middle between the two walls. The presence of the hard walls is introduced by the boundary conditions ρ(z) = 0 if z < 0 or z > L z . Equation (6) is solved in combination with equation (10) leading to density ρ(z) and order parameter S * 20 (z) profiles. The precision of this solution expressed in terms of standard deviation is 10 −5 .
The bulk order parameter S b is used both in the MF and LMF calculations. To obtain this quantity, the integral equation (5) was applied. It was also solved numerically by the iterative method, but here with a precision 10 −12 and with a step of integration over cos(θ) taken equal to 0.00125. The cut-off radii r c = 6.0/α 2 and 12.0/α 2 were used to obtain S b as well.
We compare the numerical results calculated from MF and LMF approaches with Monte-Carlo (MC) simulation results. For this purpose, a series of MC simulations in canonical ensemble [27] were performed to obtain S b , ρ(z) and S * 20 (z). A system of N p fluid particles interacting with the pair potential (1) were placed into a rectangular box of a size L x ×L y ×L z , where L x = L y = 24/α 2 and L z = 4r c if r c = 6.0/α 2 and L z = 3r c if r c = 12.0/α 2 . For L z taken in our simulations, a bulk-like region in the middle of the box is observable up to the temperature T * = 2.2. Since the bulk fluid density is set to ρ * b = 1.0, a number of fluid particles N p is set to 13824 for the case of r c = 6.0/α 2 and 20736 for the case of r c = 12.0/α 2 .
The simulations were carried out with the periodical boundary conditions applied in three dimensions in the case of a bulk fluid, and in X and Y dimensions in the case of a fluid between two hard walls.
Each simulation procedure was performed at a constant temperature and volume and consisted of three stages: 1) equilibration of a fluid in a strong field applied along Z -axis to give the fluid particles a preferential orientation, which was normal to the wall surfaces; 2) equilibration of the system obtained at the previous stage with the field switched off; 3) production of the necessary characteristics for a system obtained at the second stage with the field switched off. A criterium defining an equilibrated system is a stabilization of the total order parameter in the system, i.e., when the order parameter fluctuates around one average value during an essential number of MC steps, usually it was over 20 000 steps at least. It should be noted that in our simulations, one MC step corresponds to N p trial translational or rotational movements.

Bulk fluid
The order parameter of the bulk MS fluid is calculated at temperatures in the range of T * = 0.5 − 3.5 using the MF approximation with the different cut-off radii (figure 1). We have checked at which r * c the results tend to the case of r * c → ∞ and have tested what the effect of the cut-off radius in general is. This information is valuable to make a comparison with the computer simulations in which the cutoff radius is used. As can be seen in the case of r * c = 6.0, the MF approximation leads to the values of the order parameter lower than in the cases of r * c = 12.0 and r * c → ∞ for the whole stable nematic region up to the critical region of the nematic-isotropic phase transition, which appears at S b < 0.443 [22] ( figure 1, left-hand panel). Also, it is observed that the results for r * c = 12.0 totally coincide with the case of r * c → ∞ ( figure 1, right-hand panel). The values of the order parameter obtained from the simulations are systematically lower than in the MF approximations (figure 1, triangle symbols).  simulations. As has been shown above, this cut-off should be long enough to give results comparable to the case r * c → ∞. It is observed that all presented profiles have the same qualitative behavior, i.e., they have a distinct high maximum at the contact with the wall, then a minimum appears around z * = 0.5−1.0 and at z * > 5 one can see a convergence of ρ * (z) to ρ * b when the bulk-like region is reached. We note that the contact values obtained in the MF and LMF practically coincidethe cut off correction is negligible, while for ρ * (0) obtained from the simulations, we observe some small difference. The contact value of the obtained density profiles will be discussed more in detail later on. For this moment, we focus on ρ * (z) at distances z * > 0 and on comparison of the theoretical approaches with the computer simulations. One can see that at small distances z * , the LMF approximation agrees with the simulations better than the MF, while the MF approximation describes ρ * (z) better at distances z * around the minimum of ρ * (z) and larger. Moreover, at all temperatures considered in our study, the z * -position of the minimum of ρ * (z) is very close in the MF approximation to that in the computer simulations, while in the LMF approximation, the z * -position of the minimum is notably shifted towards the higher values.

Density and order parameter profiles
A different situation is observed for the order parameter profile S * 20 (z) of the MS fluid near a hard wall, which is presented in figure 3 as a normalized quantityS(z) = S * 20 (z)/S b . We take such a presentation because of an essential deviation of the theory from the simulations results for bulk order parame-

13608-7
ters (figure 1). Using the normalization we can consider all profiles on the same plot. Nevertheless, even with this normalization, one can observe an important quantitative and qualitative difference between not only the theory and simulations, but also between the MF and LMF approximations ( figure 3). In the LMF approximation, a minimum ofS(z) is found at z * = 0.54 and z * = 0.58 at T * = 1.3 and T * = 2.0, respectively, while in the MF approximation no minimum has been noticed at all. At the same time, the simulation results do not give any evidence of aS(z) minimum either. In this context, the MF results are similar to those obtained from simulations. However, all the considered approaches give completely different contact values ofS(z). The contact value of the order parameter obtained from the simulation is the lowest one, then an essentially largerS(z) is given by the MF approximation and the largest value of S(0) close to 1.0 is calculated from the LMF approximation. In the bulk region (far enough from the wall), all the profiles converge to 1.0 as expected. The only difference is found in the rates of this convergence, which is lower in the MF approximation than in the LMF approximation, and in the simulations it is the lowest one. Also, it is worth noting that the rate of convergence of S * 20 (z) to S b [i.e.,S(z) to 1.0] decreases with the temperature. Before discussing the contact values obtained in our study we will try to understand the behavior of the order parameter profile S * 20 (z) and how it relates to the density profile ρ * (z). First of all, we have noticed that the order parameter of fluid particles next to the wall (z * → 0) is smaller than the bulk order parameter S b . However, the density profiles (figure 2) for the same z * is higher than the bulk density.
From the knowledge of the bulk, for which the higher density leads to the higher order parameter, one can expect that the order parameter near the wall is also higher than S b . Therefore, we encounter a contradiction, since both the MF approximations and simulation predict the values of S * 20 (z) smaller than S b everywhere except the bulk-like region where S * 20 (z) is equal to S b ( figure 3). Apparently we are dealing with two competing effects: densification of fluid particles near the wall, which should increase S * 20 (z) at small z * and the absence of fluid particles beyond the wall, which should decrease S * 20 (z). Following the obtained results, one can conclude that the latter effect is more essential for the system considered in our study.

Contact values
The density and order parameter profiles of a MS fluid near a hard wall are used to calculate the corresponding contact values ρ * (0) and S * 20 (0) as functions of the temperature. These contact values can be calculated from the expressions of the contact theorem (CT) (22) and (23), which in the MF approximation reduce to the relations (24) and (25), respectively. In the LMF approximation, the contact values ρ * (0) and S * 20 (0) satisfy the relations (24) and (25) automatically due to the boundary conditions applied in the solution of the corresponding differential equation. In the MF approximation, the situation is different for ρ * (0) and S * 20 (0). As it was shown in [23], the contact value for the density profile ρ * (0) satisfies the relation (24), although it has not been proven that the contact value of the order parameter S * 20 (0) obtained in the MF approximation should satisfy the relation (25). First we consider the temperature dependencies of the contact value of density profiles obtained with the different cut-off radii r * c = 6.0 and 12.0 (figure 4). As can be seen, the CT leads to the same result as the MF approximation. A small difference between the LMF and the MF approximation (or the CT) appears when the cut-off radius is equal to r * c = 6.0 ( figure 4, left-hand panel). If the cut-off radius is increased to r * c = 12.0, equivalent results are obtained in all of the approximations (figure 4, righthand panel). At the same time, the agreement between the theoretical approaches and the simulations is mostly qualitative. We would like to draw attention to the non-monotonous behavior of the contact value ρ * (0) with a distinct minimum observed both in the simulations and in the theoretical predictions. To understand this interesting effect, we should analyze the temperature dependence of the density contact value more in detail.
As it is seen in figure 4, at low temperatures, ρ * (0) is large and lowers as the temperature increases.
This is related to the reduction of the repulsive contribution of the pair potential (1). Since in our model a soft repulsive interaction is used, it becomes weaker if the temperature increases. It should be noted that the attractive part of the pair potential in our model is rather small and mainly affects orientational properties of fluid particles. At the same time, a relative orientation can indirectly strengthen the fluidfluid repulsion contribution by reducing the attractive interaction term. Therefore, at some point in the solely, models with a pair potential consisting of a soft-core term combined with a Maier-Saupe attractive potential. In the case of a hard-core type of repulsive interaction, the temperature dependence can be opposite and the effect of orientational ordering can be completely different.
While the contact value of the density profile is rather understandable, the behavior of the contact value of the order parameter profile is not so clear. First of all, as it has been already shown, there is an essential inconsistency between the MF and LMF approximations ( figure 3). There is also a problem if one compares the contact values S * 20 (0) obtained in the MF and LMF approximations and those calculated from the CT theorem. As can be seen in figure 5, the MF results significantly differ from those of the LMF and the CT. A perfect agreement of the LMF and the CT appears due to the definition of the boundary conditions used in the LMF approximation, which are taken to fit the CT theorem (25).

13608-9
r * c = 12.0, the contact values of the order parameter of a MS fluid obtained from the CT and the LMF totally coincide. At the same time, the MF approximation gives much lower S * 20 (0), although it is closer to the simulation results than the CT and LMF. We assume that the difference between the CT and the MF approaches is related to a poor description of the bulk partial pressure P (Ω 1n ). In order to improve the results obtained from the CT (or the LMF approximation), one should take into account higher order terms, for instance the Gaussian fluctuations [22]. Also, it should be noted that in contrast to the CT and LMF approximation, S * 20 (0) obtained from the MF approximation and the simulations decays to zero before the critical region is reached. Thus, the convex part of S * 20 (0) is absent for them. It means that in the critical region there is no orientational ordering of fluid particles at a contact with the wall.

Conclusions
Within the framework of the density field theory, we have formulated the mean field (MF) approximation as a starting point for the theoretical description of a nematic fluid at a hard wall. Using the developed approach we have investigated the density and order parameter profiles of a confined Maier-Saupe nematogenic fluid with an isotropic Yukawa-like repulsion. Theoretical predictions have been compared with analytical results obtained within the framework of the linearized mean field (LMF) approximation [23]. For the density profile, the results obtained in the MF and LMF approximations are in qualitative agreement with computer simulations data. For the order parameter profile, the results of the MF and LMF approximations have qualitatively different behaviors. In the LMF approximation, a minimum of the order parameter profile is present, while in the MF approximation no minimum has been observed. The results of computer simulations do not give any evidence of the existence of a minimum in the z-dependence of the order parameter profile either. In this context, the MF and computer simulations results are qualitatively similar. This contradiction in the description of the density and order parameter profiles within the framework of the LMF approximation is connected with the problem of the respective contact theorems in forms (24) and (25) which are used as boundary conditions in the solution of the corresponding system of differential equations. As it was shown in reference [23], the contact value of the density profile ρ(0) satisfies the relation (24) but the validity of the relation (25) for the contact value of the order parameter S 20 (0) is not evident. Moreover, the comparison with the results of computer simulations in figure 5 shows that relation (25) is probably incorrect. As we can see from figure 5, the LMF approximation can give a better result if we change the boundary condition for the order parameter profile by correcting the contact value for S(0) obtained within the framework of the MF approximation. We hope that the considered problem can be better understood by going beyond the MF approximation and including a contribution from fluctuations.
The temperature dependencies of the contact values of the density and order parameter profiles have been analyzed more in detail. We have found a non-monotonous behavior of the density contact value ρ(0) as a function of the temperature with a distinct minimum observed in both theoretical predictions and computer simulations. This non-monotonous temperature dependence of ρ(0) is explained by the competition of the soft isotropic repulsive and the soft anisotropic attractive contributions. For the contact value of the order parameter, we have observed a monotonous decrease with an increase of temperature. Both the simulations and the MF results indicate that there is no orientational ordering of fluid particles at the contact with the wall when the fluid is in the critical region.
The results presented in this paper have been obtained within the framework of the MF approximation. We note that the agreement between theoretical predictions and computer simulation data is mostly qualitative. For a better theoretical description one should take fluctuations into account. For the bulk properties of the model under consideration, the influence of the contribution from fluctuations was already discussed in reference [22]. It was shown that the singlet distribution function reduces to the form (5) with a change of κ 2 2 to κ 2 2 t , where t = 1 − βA 0 α 2 2 /(λ 0 + λ 2 ). It was shown that the temperature trend of deviation between the order parameter S b calculated with fluctuations included and the MF value is the same as that between computer simulations and the MF value presented in figure 1. In the next paper we plan to include the contribution from fluctuations in the description of a MS nematogenic fluid at a hard wall similar to the way it was done for isotropic Yukawa fluids [18,19].