Structure formation of colloids in nematic liquid crystals

We investigated the behaviour of colloidal particles suspended in nematic liquid crystals. These colloidal particles interact through elastic deformation of the nematic director field which can result in nontrivial collective behavior, leading to the formation of spatially modulated structures. In this paper, the formation of lattice structures is described both by computer simulations and by analytical theory. Effective interactions of the pairs of spherical macroparticles suspended in nematic liquid crystals have been suggested by many authors. Using these pairwise interactions, spatial structures are obtained by means of dynamic simulations. We have suggested a number of possible structures, which may be formed in multi-macroparticle systems. Regions of temperatures and concentrations are determined in which such a structure might appear.


Introduction
Colloidal structures formed by hard particles suspended in liquids represent a widely abundant state of matter and have attracted much interest in science, technology and medicine [1,2].Among those colloidal systems, a liquid crystal mesophase containing macroscopic inclusions of foreign substance is a particularly interesting example.The physical properties of such systems are thoroughly studied in [2][3][4][5].However, as pointed out in these studies, there still remain many interesting problems to be solved.Liquid crystals (LC) are anisotropic liquids macroscopically characterized by the field of the director which points in the local average direction of the constituting molecules.When macroscopic particles are introduced into LC's, the electro-optical properties of the system are considerably effected, which is often due to the new spatial structure formed by the particles.
To describe the conditions of formation and the properties of these macroparticle structures, one should take into account both aspects of the particle interaction: those associated with physical properties of macroparticles (such as the electrical charge) and those caused by the deformation of the director elastic field.Interactions among the dispersed particles resulting from the deformation of the nematic director field were investigated by many authors [3][4][5][6][7][8].The physical mechanism is that the dispersed macroparticles cause a distortion of the director distribution and via this deformation of the elastic field other macroparticles are affected.This gives rise to an effective pair-wise interaction between two macroscopic particles dispersed in nematic LC(NLC).In other words, the interaction between a pair of macroparticles in NLC can be effectively described without directly treating the background director field.It is interesting to note that there are two different treatments concerning the surface anchoring.In one approach, the director field near the particle is determined solely by the director distribution at the particle surface and does not depend on the concentration of other particles in NLC.In this case, the strong deformation of the director field yields defects [6,7].In the other approach, the director distribution at the surface of an individual particle is determined consistently with the director field deformed by other surrounding particles [8].In the above two cases, the effective pair interactions of sphere particles are of quadrupole type: φ 2 (r ij ) ∝ r −5 ij .However, the dependence on the angle θ ij of r ij against the far field director n is different for [6] and [7,8].A completely different approach which avoids treating the anchoring at the surface of the sphere particles [9] also leads to a quadrupole interaction with the same angle dependence as [7,8].Thus, we can say that the quadrupole type interaction exists in general as an effective pair potential among the sphere particles dispersed in NLC.As a matter of fact, it has been shown that a quadrupole interaction always exists for sphere particles no matter whether the anchoring is normal or planar or whether the defects exist [10].In addition to the quadrupole type interaction, dipole type interaction emerges when the mirror symmetry in the distribution of the director around the particle is broken [10].Since a quadrupole interaction always exists for the spheres embedded in NLC, the states of multi-body systems obtained by using such effective pair potentials are important as reference states.Depending on the anchoring strength, higher body interactions might exist such as three body φ 3 (r i , r j , r k ) or, more generally, n-body interactions φ n (r i , r j , . . ., r n ), which, in essence, strongly depend on the specific configuration of the embedded particles.No general expression exists for such higher body interactions.However, when the anchoring is weak we can estimate the leading terms as φ 3 (r 12 , r 23 , r 31 ) ∝ α 2 γ, where anchoring and surface energies are respectively proportional to α = W r 0 /K and γ = W r 2 0 .Since φ 2 (r ij ) ∝ αγ, we can expect that the pair interaction φ 2 becomes dominant when α 1.In such situations, the structures observed in this work might actually appear.
In this paper, we treat the quadrupole interaction of spherical particles obtained in papers [6][7][8][9].The only difference with φ 2 in [6] and [7][8][9] is ψ(θ ij ); the dependence on the angle between r ij and n where the interaction of spherical particles i and j introduced into NLC may be either attractive or repulsive depending on the angle θ ij .Therefore, a question naturally arises, if it is possible to form a new spatial-periodic structure in such a system of interacting particles, and what the conditions are of such a structure if the answer is positive.The purpose of this paper is to give both computer-simulated and analytical descriptions of the conditions and the occurrence of such structures.At first we will carry out computer simulations of the systems of multi-particles cooperating through pairwise interactions.Possible structures will be obtained through dynamically finding the stable position of the particles.There are several cases when there is an obvious structure formation.These structures are formed as a result of the collective behavior which is described by pairwise additive interactions.It will be shown that the collective behavior of the particles depends on ψ(θ ij ).Next, we will give the analytical solution of the problem.We will present a possible approach to the description of the structure formation which has been suggested in papers [6][7][8].External physical conditions (concentration, temperature) and geometrical characteristics will be determined concerning the realization of the formation of spatial-periodic structures.

Computer simulation and results
In the following section, we describe the model and simulation method and discuss our results.We consider a system of spherical foreign macroparticles in a nematic mesophase.The pairwise interaction of the particles can be expressed as a function of distance, angle and anchoring strength, implicitly including the effect of collective deformation of the background elastic director field.To discuss the qualitative difference of the effective potentials we investigate only their angle dependence.The potentials we use are of LT [7][8][9] form: and of RT [6] form: where r ij is the distance between particles i and j, r 0 is the radius of the particle, and θ is the angle of r ij against the far-field nematic director.We take the director to be parallel to the z-direction; cos θ = ( r ij • z ) where z is the unit vector in z-direction.
The first term (r 0 /r ij ) −12 in equations ( 2) and ( 3) is added to avoid the collapsing of the particles into each other when interparticle forces are attractive, and represents the excluded volume effect of the colloidal particles.We set the mass, r 0 and to unity, respectively.For the simulation results reported here, we use initial configuration of a stacked triangular lattice structure in nearly a rectangular parallelepiped volume with various lattice spacings.No boundary conditions have been applied.We choose to let the system relax and determine its density spontaneously to investigate purely the effect of the potentials.The system was relaxed through particle dynamics simulations using 6-value 2nd order Gear predictor-corrector method [14].
We start with the same initial conditions for both interactions ( 2) and ( 3), for the same particle numbers to see the different nature of the interactions.We used two methods of simulation and confirm that the structures obtained in this paper are stable.In both methods, we use random initial velocities which obey the Maxwell distribution (corresponding to the Brownian motion of the particles).In the first method, we relax the system by letting the particles obey the Newtonian equations of motions.The system will spontaneously evolve to a state where the value of potential decreases and the position of particles stabilizes.Moreover, in the second method, we scale the particle velocities to maintain the kinetic energy of the system.In such a case, the temperature of the system appear to be the fixed kinetic energy of the particles.In these cases, the thermal fluctuations of the particles corresponding to the temperature of the system go on.However, the system evolves to a state where the force acting on the particles becomes small enough (of the order of 10 −4 ) such that we can assume that the system has reached a stable state.From this second method, we know that a critical temperature exists having to deal with the structural formation.Simulations of temperature T = 1[ ] = A/r 5 0 result in stable structures (reported in this paper), while the temperature increase to T = 2 results in the particles dispersed in a vast volume.In both cases, the particles disperse randomly in a vast space when no stable configuration could be found.
The results of computer simulations for 10 particles in the system are presented in figures 1 and 2 for both potentials (LT model and RT model).We would first call your attention to the structure formed by these macroparticles in the system.Chain-like structures, which are similar to the experimentally observed ones [4,11], can be formed when the number of particles is small (figures 1 and 2).The minimum of the pair's potential energy gives an angle of cos 2 θ = 3/7 in case of equation ( 2) and cos 2 θ = 4/7 in case of equation (3).So, both LT and RT models cannot give rise to an isolated single chain parallel to the director of colloidal particles since the pairwise interactions only yield a repulsive force when two particles i and j align along the director (θ = 0).However, the clusters can give rise to chain structures due to the neighboring particles.The particles in oblique direction in the neighboring chains are the source of attractive forces which stabilize the chain structure in this case.As evident in figure 1(b) and figure 5(b), the chains are approximately parallel to the director.In experiments, chain-like structures are observed parallel to the director for N D phase and oblique for N c phase [4].In these experiments, the chain-like structures of colloids appearing parallel to the director seem to be fairly isolated from each other.Quadrupole-quadrupole interactions that we treat here cannot give isolated chains parallel to the director, but they are possible to form clusters including chainlike structures.
The pairwise interaction energy given by the RT model gives rise to a periodic structure when there is a small number of particles in the system.However, when    the number of particles in the system increase, the particles distribute much more randomly in space.In contrast, the LT model gives rise to rather well-defined layer structures even for larger numbers of particles (figure 3).It seems that the LT model represents the structural formation more satisfactorily in many-body systems, minding that the initial configuration and the conditions of simulations are identical with the same number of particles.Since the potential (2) should satisfactorily describe the collective behavior of the colloids, the structure formation itself is not a surprise.We would like to call your attention to the hexagonal distribution of the particles in the layers perpendicular to the director shown in figure 3(a).Similar droplet-like structure with the layers of hexagonal structure can be seen with the other cluster sizes.In these layer structures, the periodic modulation in parallel and in perpendicular directions gives the ratio R /R ⊥ = 3/4 (or cos θ = 3/7).This ratio is where ∂U 1 /∂R = ∂U 1 /∂R ⊥ = 0.In figure 4, we show the pair distribution function in the directions perpendicular (a) and parallel (b) to the director with solid lines.Strong modulation in the direction parallel to the director in figure 4(b) shows that there is a layered structure.In the direction perpendicular to the director (figure 4(a)), periodic modulation can be discerned for several particle neighbours.Both pair distribution functions decay towards the cluster size.In contrast, the RT potential does not form a well defined cluster structure as shown in figure 5. Note the difference in scale for figures 4 and 5.The RT potential leads to a random configuration dispersed more widely in space (more than three times in length compared to the well defined cluster structure obtained in the LT model).Also, note that the cluster obtained from the LT potential is more droplet-like (the length in the direction parallel and perpendicular does not differ greatly), while the RT potential leads to the particles dispersed more widely in the direction perpendicular to the director.

Theoretical model
In the following section, we give an analytical description which shows that the structure formation observed in simulation can occur for certain concentrations and temperatures.We can employ an approach similar to molecular field approximations [8].From a physical point of view it is obvious, for reasons analogous to those set forth in [13], that the macroparticle distribution obtained from the interactions (2) and ( 3) can be inhomogeneous.The results in paper [8] illustrate the possibility of describing the system of interacting particles by means of a statistical theory.The proposed approach provides a unified statistical description of the systems of interacting particles with spatially inhomogeneous distribution.It also allows one to consider the collective behavior of the particles with attractive and repulsive interactions, and determine the conditions of the cluster formation in the model system.
In our analytical description, we consider the system of foreign macroparticles suspended in the mesophase as the Ising model lattice gas without fixing the basis of the initial arbitrary lattice.This treatment is completely analogous to the model of substitution of interstitial solid solution [13] and we will employ the approach proposed there.The free energy of solution of macroparticles in a liquid crystal in a self-consistent field under approximation of binary interaction, may be written in the form: where U ( r − r ) is the interaction potential, f ( r ) is the probability that the macroparticles fill the initial arbitrary lattice, µ is the chemical potential, and T is the absolute temperature.Since the particles cannot occupy the same space position, it is natural that equation ( 4) coincides with the entropy of Fermi particles.The minimization of equation ( 4) yields the self-consistent field solution for f ( r ) i.e.
Each set of thermodynamic functions of state correspond to a solution that describes some phase of droplet arrangement.If their distribution can be inhomogeneous, then the solution (5) serves as the stable phase associated with the interaction character and temperature.An efficient way to obtain such a distribution was proposed in [13].If the solution of macroparticles is disordered, then by definition (the mean value) of f ( r ) = c, where c is the relative macroparticle concentration.The concentration inhomogeneity gives rise to an additional term f ( r ) = c ± ϕ( r ) where ϕ ( r ) is the change in the probability distribution function of the macroparticles.
If the concentration inhomogeneities are smooth and their scale is longer than the interparticle distance, the quantity may be interpreted as the change in macroparticle composition.When passing to a continuum description, we can write the free energy increment, associated with the inhomogeneous macroparticle distribution in terms of the power series expansion by using the long-wavelength expansion of the concentration i.e.
where ρ = r − r .In this case we may rewrite the change in free energy in the form: where ρ and ρ ⊥ are the components of the intermacroparticle distance vector along and perpendicular to the director.Thus we see that the minimum of this function realizes a spatially inhomogeneous macroparticle distribution provided that the sign satisfies some relation and the values of coefficients determined by the intermacroparticle interaction potential are given by equation (8).In order to reveal the conditions under which the homogeneous macroparticle distribution becomes unstable, we have to calculate all the coefficients.This can always be done; however, for the sake of simplicity and clearness we consider a cylindrical sample with the base radius L and height d, the boundary condition being that the director is normal to the cylinder base.Then, for the interaction energy given by equation ( 2) we find all the coefficients in the cylindrical coordinate system [8].In order to reveal the conditions under which the homogeneous macroparticle distribution becomes unstable, and the temperature of the phase transition to new states may be determined from the following relation: The critical temperature of loss of stability of homogeneous macroparticle distribution is determined and described by the following expression where A = r 5 0 is the coefficient in the energy of interaction between foreign spherical macroparticles in nematic liquid crystals, and can be related to the Frank elastic constant K and anchoring energy W by A = (6/πK)(4πW r 4 0 /15) 2 .The first term in equation ( 2) gives only a negligible contribution to equation (10) and will be omitted in the discussion that follows.Formation of such structures may be treated as a first-order phase transition since the cluster formation can be interpreted as the formation of a new phase.The approach makes it possible to describe a system of interacting particles using a statistical method by taking into account the spatially-inhomogeneous distribution of particles [13].In general, we can describe the possibility of formation of 3-dimensional structures and other types of periodical modulation in macroparticle distributions.Here we will consider the case where concentration can be described by a wave form, and assume that the uniform lattice configuration needed to form this concentration wave is formed in accordance with the new equation.We write the concentration wave in the form of Fourier transform ϕ(k) = drϕ(r) exp(ikr) (11) of an arbitrary function ϕ(r) defined in a finite volume with periodic boundary conditions.In this case the relation which determines the inhomogeneity of macroparticle distribution is as follows: Thus, we judge that the homogeneous distribution of macroparticles dispersed in nematic liquid crystals can be unstable in a limited temperature region in which the mesophase exists, and spatially modulated distribution can be generated with wave vectors k and k ⊥ .The criterion of instability given by the relation equation ( 12) can be interpreted as a condition for the formation of spatially non-uniform distribution at a given temperature, which depends on the concentration of macroparticles and the wave vector of the new structure.The coefficients ξ , ξ ⊥ calculated for cylindrical boundary condition of radius L and height d will lead to When T = T c , the value k /k ⊥ will diverge under the parity L 2 = 2d 2 .Using this condition L 2 = 2d 2 , equation (10) will become as follows: where N is the number of the particles in volume V = L 2 d.
Using the above analysis, comparison with simulation results is possible.Figures (3) and ( 4) approximately fulfill the relation L 2 = 2d 2 .So, we can use equation (14) and the values obtained from simulation.From figure 4(b), the layer distance is d = 5.2 r 0 in the system of N = 128.By using these values we get r 0 = 1.3 µm from equation ( 14) , if we assume that structural transition occurs at T c = 300 K and the anchoring coefficient and elastic constant are W = 10 −6 J/m 2 and K = 10 −11 N, respectively.

Conclusions
We presented computer simulated and analytical results of possible structures for suspension of macroparticles dispersed in liquid crystals.We place special emphasis on cluster structures which are formed in the systems of spherical macroparticles with perpendicular director orientation on the surface by using pair potentials ( 2) and (3).We have shown that lattice structure formation can occur in liquid crystals for realistic values of temperature and macroparticle concentration.
As it is known, a system of spherical particles having spherical-symmetric potentials of interaction forms crystal structures whose pattern depends on the range of the interaction potential, the temperature, and the particles' concentration [15].A wide variety of crystal structures is possible if the forces of interaction possess anisotropy.Such an interaction can also be realized in the system of spherical macroparticles, cooperating via deformation of the elastic field of the director in a liquid crystal.Another important property of such systems is that they are highly sensitive to weak external effects.It enables us to exert action upon the formation of the structures and their transmutations.Therefore, the macroparticles introduced into liquid crystals represent a system for research of conditions of formation of possible crystalline structures.In addition, this system is advantageous in that they are visually observable and one can directly observe all changes of the structures.Thus, we have a perfect modelling system for determining the physical parameters of the medium, necessary for effectively controlling the structural formation.For this reason, colloid science and other fields of physics sometimes refer to a soft condensed matter.The effect of external fields is, therefore, relatively strong and can be used to manipulate the crystallization process.Quantitative real-size measurements of these structures are possible via the employment of high resolution and sectioning power of confocal microscopy.Here we consider a system whose interactions are presumably wellunderstood and in which the phase transition can be studied in the framework of macroparticle resolution i.e., the model system permitting the study of microscopic mechanisms of the structural phase transitions.All we assume is that the pair potential is appropriate in describing the interaction in a multi-body system.Under this assumption, the experimentally determined phase diagram for bulk colloidal suspension has been found to resemble that obtained from dynamic simulation of particles interacting via elastic deformation distribution of the director field.Quantitative discrepancies have been quite reasonably ascribed to the uncertainty in the experimentally controlled parameters.Even the elastic properties of colloidal crystal seem to be consistent with the theory.However, this qualitative agreement should not be the reason for complacency.The still evolving theory of phase transitions tells us that, qualitatively, the phase behavior of materials depends to a degree on the shape of the pair potential.In this paper, we treat colloidal systems in the mesophase by using pair-wise potential to effectively describe the elastic deformation of nematics.When the pair interaction potentials are of a shape suggested in the papers [8,6,7], we describe the possibility of the existence of a supermolecular structure, which may be formed in the system of foreign macroparticles in a liquid crystal.Numerical results are consistent with a liquid-structure theory recently proposed to address the phase separation in an interacting colloidal suspension.Consideration of multi-body effects or other interaction, for example dipole-dipole interaction, remains to be explored and should be an important step towards the complete theory of colloidal interactions.
¢¡ ¤£ ¦¥ ¨ § © £ !§ " # # § £ ¦$ %" &¡ (' )¡ 10 32 54 6© 7© 7 8 ¦¥ 9 $ @ $ !A £ ¦4 2 5" B$ !A C" B£ ¦$ B ' D BA E GF IH 6F P 5Q (R TS VU W )X `Y ¦F Ia bF c )d e gf hS X `i 5F Ia bF p gd q sr ut e W )X `v wF Ix yd e d !q s S 3U S 1 1 ¤ 1 w & g B 3 ed u ` gf d u gh i g h 3j uk uk 1 ml in uj u & o p 1q i r W sf q t u ev eu w hu w G x v g hv sy `z y T{ | hh i } ~q iv ` k #k gl ul s 1v } ~ u p uh u gt ip 1 gy B e hv ` #n %y B e w 1 g ` 6 | gt 1 q iv g G q t u ev eu w hu ¦q e gq u g q g 1 3 v G | gt 1 #¡ e¢ u £p 1 g g q i 3 w g i} q iv g G¤ h it u ev gq 1 D sq g !eu ev i¤ eq iv g 5h # hv g hh #v it u h e u g g w ( ¥ w g i} q i ¦¤ h it u ev gq i hv ¤ x g i¢ u § (u e|¤ gh u x p uh ©¡ ¦q | w D 1h # !# ª § «p ur q g !eu ev i¤ eq 1 ¥ 5 v gh u g 3u 1h i ` #¬ G 9 ¡ t ip uh #v i¤ ev gq i § ¤ eu ev q eu eh #v g e 1gq ) w g g 3u ev g q ¤p 1 g g q i 3 g hp uh #v g g &r ¤ ev | Dw D 1h # gq iq er )p uh u gt u 1h u 1 g ` 1 h 1ḧ 1 gq iv g t u eh 3u w h & ª t u h eu wu e gw D 1h # gq iq er ® 3h i eu e w 1 v g t u eh u w h 1p uv it e¢ gu eḧt er Gx g ) | 1p 1 1 ¥ 1 § ¯r u ° w 1 p 1± § (u gh #q 1 ¥ T | e 1 § ¤ gq iq er 1 u h g ° gq e e 1 u ev i¤ eq i} wu i 1h # } v g w 1h #v it u 1 © § °¤ ev ¦h i gq i ² D x g gp uh u 1p 1q 1 gq i w D g 3u ev g q i p gh #q i ¥ x g i¢ u } © !g hh u g¤ h it u ev gq 1 | q g !eu ev i¤ eq iv g °h # hv g ( hh #v it u h e u g g V eu ` | 1 v gq e g ¤ eq 1 ¥ T | e 1 § ¤ gq iq er y u eh #v g !gq p uh u gt u 1h u 1 3t u eh 3u w h #v ` ³ | # p uh u 1p 1q 1 gq yh ir ¨ ¡ 1v g v g °t u eh 3u w h & Vr u h 3 ¡ hu e!w D 1h # eu ev ) h 1¨u ev # !g hh u g¤ h it u ev gq i w 1 Vt iv it u g v ` v x q e i¤ gq gd u it u e ´u g p gh i eu w h s w w 1q ª gq u h i ª ` g r u hv g (u h g h 3t u eh 3u w h #v y ¡ hu e v gq iv g eu ev ` µ ` ¶ • si¹ |º g» h¼ i ¶ |¹ |º u½ B¾ ¿ À gÁ &À #Â Ã Ä #Å AE #Ç È É iÊ ´Ä &Ë AE uÌ ´Í &Ì Ã ¿ ©Ì ¿ VÍ &Ä 1Î iÊ ´É gÁ 1Ä Ï gÐ (Ñ Ò ¾ eÓ Ô `Õ Ö × uÕ 6Ø IÙ ©Å gÚ Û VÕ ÝÜ 1× uÕ 6Ø mÞ ¥Å uÓ Ô `Õ Ö × uÕ Ýß (à

Figure 1 .
Figure 1.Structure obtained for LT potential with N = 10 (a) Director is perpendicular to the plane of the figure (bird view) (b) Director is parallel to the vertical axis (side view).

Figure 2 .Figure 3 .
Figure 2. Structure obtained for RT potential with N = 10 (a) bird view (b) side view.

Figure 4 .
Figure 4. Pair distribution function (2D) in the direction (a) perpendicular and (b) parallel to the director obtained for LT potential with N = 128.In addition, the ordinary 3D pair distribution function is shown in a broken curve in (a).

Figure 5 .
Figure 5. Pair distribution function (2D) in the direction (a) perpendicular and (b) parallel to the director obtained for RT potential with N = 128.In addition, the ordinary 3D pair distribution function is shown in a broken curve in (a).
` D{ | x e 1q | g B 3 ed u Dj uk u 1 m i un u & o p 1q i r ) u | & ¤ & g ¤ 3 ` & &