Field Induced Memory Effects in Random Nematics

We studied numerically external field induced memory effects in randomly perturbed nematic liquid crystals. Random anisotropy nematic-type lattice model was used. The impurities imposing orientational disorder were randomly spatially distributed with the concentration p below the percolation threshold. Simulations were carried for finite temperatures, where we varied p, interaction strength between LCmolecules, and impurities and external fieldB. In the {B, T} planewe determined lines separating short range— quasi long range and quasi long range—long range order. Furthermore, crossover regime separating external field and random field dominated regime was estimated. We calculated remanent nematic ordering in samples at B = 0 as a function of the previously experienced external field strength B.


Introduction
Impact of weak disorder on statistical properties in phases exhibiting continuous symmetry breaking is of interest for different branches in physics [1][2][3].Such systems are of particular interest due to memory effects which they can exhibit [2,4,5].Pioneering studies have been done in magnetic systems [1,6].Lately it has been shown that various liquid crystal (LC) phases are extremely adequate for such studies [7].Namely, LCs are relatively easily experimentally accessible, they display rich diversity of different configurations, and there exists plenty of different means how qualitatively different types of disorder could be imposed to LCs.Consequently, LCs are ideal testing ground to investigate fundamental problems related to impact of weak disorder on phases reached via a symmetry breaking phase transition.In addition, randomly perturbed LCs are also technologically relevant, specially due to their unique combination of liquid character, softness, and optical anisotropy.Randomly perturbed LCs could be exploited for various devices, particularly related to electrooptic applications.
Uniaxial nematic phase represents the simplest LC configuration exhibiting only orientational ordering [8].At mesoscopic level this ordering is commonly described by the nematic director field ⃗  and the scalar nematic order parameter  [8].The unit pseudovector ⃗  points along local average uniaxial LC orientational ordering, exhibiting headto-tail symmetry (i.e., states ± ⃗  are physically equivalent).On the other hand  measures extent of fluctuations about ⃗ .Absence of thermal fluctuations in local orientational order are fingerprinted in  = 1.In bulk equilibrium nematic phase is spatially homogeneous where ⃗  is aligned along a single symmetry direction.In thermotropic LCs, to which we restrict our study to, nematic long range ordering (LRO) is obtained via the first-order phase transition from the isotropic phase, by decreasing temperature.The latter phase corresponds to ordinary liquid ordering, where  = 0.
There are several different ways how random-type disorder could be experimentally imposed to LCs [7].For such purposes LCs are either immersed in various porous materials (e.g., aerogels, Vycor glasses, Russian glasses, controlled pore glasses, and millipore membranes. ..) or by mixing them with different nanoparticles.In the latter case aerosil nanoparticles are commonly used [9].Namely, by changing their concentration, three qualitatively different types of disorders could be realized [10].
Resulting randomly perturbed nematic structures depend on disorder strength.According to the Imry-Ma Advances in Condensed Matter Physics theorem [1,6], one expects that weak random field-type disorder stabilizes a domain-type pattern characterized by a single length   obeying the scaling law   ∼  −2/(4−) .Here  stands for the disorder strength and  is the space dimensionality.The resulting phase should possess short range order (SRO).However, several subsequent studies suggest that the quasi long range order (QLRO) or even LRO might be established [2] instead.
In general, randomly perturbed nematic structures possessing SRO or QLRO exhibit glassy features [2].Namely, strong enough disorder could stabilize enormous large variety of different metastable structures in free energy landscape.Recent experimental and simulation studies in different randomly perturbed nematics support this expectation [4,9,[11][12][13].For example, in LC-aerosil mixtures, it was demonstrated that LCs could memorize external electric or magnetic field induced orientation [4,5,14,15].Note that systems exposed to relatively strong disorder are known to exhibit several universalities, independent of systems' microscopic details [2].However, in case of relatively weak disorder system, specific details are expected to emerge, which further increases diversity of potentially possible phenomena.In addition, glassy features in nematics could be apparently different with respect to magnetic analogues due to nematic head-to-tail symmetry, which is not observed in magnetism.Consequently, perturbed nematics structures in general display complex tangles of topologically stable line disclinations characterized by winding number  = −1/2 [5].On the contrary, line defects could not be topologically stabilized in magnetic systems.Namely, local magnetization could exhibit only line defects with winding number  = 1, which could avoid singularity in orientational ordering by escaping in "third dimension." Furthermore, in LCs, thermal fluctuations play relatively important role due to their soft character.Consequently, energy barriers in free energy landscape could be relatively overcrossed easier with respect to magnetic analogues.
In the present paper we study disorder enabled memory effects in randomly perturbed LCs.Of our particular interest is external field driven memorization.Our study is numerical, where we use lattice version of modified Random Anisotropy Nematic model [16][17][18].The plan of the paper is as follows.
In Section 2 we present the model.In Section 3 results of our simulations are presented.In the last section we summarize results.Some technical details related to numerical method and mesoscopic interpretation of results are collected in two appendices.

Model and Numerical Approach
We use lattice-type simulations where the cubic simulation box of volume  3 contains  = (/ 0 ) 3 sites.Lengths  and  0 represent the characteristic linear size of a system and the shortest distance between adjacent lattice sites, respectively.Typically, we set / 0 = 70.Nematic LC orientational order at each lattice site is presented by the unit vector ⃗   , to which we henceforth refer to as a nematic spin.The so-called head-to-tail symmetry ± ⃗   is considered to mimic typical nematic LC behavior at the mesoscopic level.In addition, we suppose that each site might be populated by a rod-like object with probability .Their local quenched orientations ⃗   are determined randomly, exhibiting isotropic orientational distribution, where the head-to-tail symmetry and | ⃗   | = 1 are imposed.We refer to these pseudovectors as impurity spins.Note that sites occupied by impurities are not accessible to nematic spins.
Using random number generator we first set static impurity positions and their orientations.For such initial configurations we then determine fixed point configurations of nematic spins, corresponding to metastable or globally stable structures.
2.1.Hamiltonian.The interactions within a system are described by the following Hamiltonian [5,11,19,20]: ( The first part of (1) describes short range interaction among nematic spins, where   =  > 0 for neighboring spins (each site has 6 first neighbors) and otherwise   = 0.The factor 1/2 is added in order to avoid double counting of pair-interaction contributions.The second term takes into account presence of a global uniform field ⃗  =  ⃗   , which we align along the -axis of the Cartesian coordinate system {, , }.The last term takes into account short range interaction between an impurity at a site  and a nematic spin at a site  of strength  > 0. Therefore, if impurity is not present at a site  then   = 0.For such choice of interactions a nematic spin tends to be oriented parallel to a neighboring nematic or impurity spin.Furthermore, a finite external field tends to align nematic spins parallel to its direction.A resulting configuration compromises couplings represented by interaction strengths , , and .Simulation details are described in Appendix A.
In simulations we use three qualitatively different histories of systems.We imposed either temperaturequenched history (TQH), field-quenched history (FQH), or temperature-annealed history (AH).The first two of these limits represent extreme histories encountered in typical experimental studies.In case of TQH (FQH) we quench system from isotropic initial distribution (homogeneously alignment along ⃗ ) of nematic spins.In case of AH we calculate a fixed point configuration { ⃗   } at a given temperature originating from an initial fixed point nematic structure which was calculated at a slightly different temperature.Simulations were repeated  rep ∼ 10 times for a given set of parameters (i.e., , , and a chosen initial condition).

Monitored Quantities.
From obtained fixed point configurations we calculate the global order parameter  () and local order parameter tensor  () .They estimate a global and local degree of orientational ordering, respectively.We define them as follows: where the index  runs over the lattice sites occupied by nematic spins,  runs over all the neighbours (their number is   ) of a site ,  is the identity tensor, and ⊗ marks the tensor product.The average local scalar order parameter  is defined as the largest eigenvector of  () .On the other hand the average global scalar order parameter  is defined as the largest eigenvector of  () .We further define the orientational two-point correlation function as [11]: where ⟨⋅ ⋅ ⋅ ⟩ is the statistical average of the scalar products of nematic spin pairs separated by a distance In case of short range order (SRO) or quasi long range order (QLRO), it holds ( → ∞) → 0. In case of long range order (LRO), it follows ( → ∞) →  ∞ ∼  2 .In general, one expects an exponential decay towards a saturated value of () on increasing  for both LRO and SRO.On the other hand, for QLRO, algebraic decay of correlations is expected; that is, () ∝  − .
In order to obtain structural details from () for a finite system the correlation function is fitted using the empirical ansatz [11]: where , , , and  are adjustable parameters.The quantity  estimates an average linear size  = 1/ of relatively wellcorrelated regions, referred to as domains.Note that, in case of QLRO, the decay of correlations with distance is relatively weak, and finite-size effects are expected to be important.Consequently, the power law coefficient  determining algebraic decay of correlations is estimated by the following equation [11,21,22]: where  ∼ /2, see (4).In case of SRO and LRO, one expects  ∼ 3/2 and  ∼ 0, respectively, and QLRO is fingerprinted by 0 <  < 3/2 [11,23].

Results
We calculated fixed point configurations of nematic spins on varying parameters {, , , } and histories of systems in three-dimensional systems.We limited to concentrations of impurities below the percolation threshold, corresponding to the critical concentration   ∼ 0.3 [14].Of particular interest is impact of external field  on degree of ordering and memory effects.
In Figure 1 we plot a typical temperature evolution of average degree of ordering () in a temperature loop using AH absence of external field.We started simulations at a relatively high temperature  =   >  IN , where nematic spins exhibit essentially isotropic ordering.Here  IN marks the isotropic-nematic phase transition temperature of a bulk sample.In our simulations we obtain  IN ∼ 1.21 which is in line with existing simulation results using similar modelling [11].Then we were gradually decreasing temperature until the temperature   = 0.5 was reached, corresponding to conditions deep in the nematic phase.Afterwards, we were gradually increasing temperature till  =   .The plot () reveals that, for example, for {, } = {0.25,0.1}, the system displays thermal hysteresis, suggesting first-order phase transition.Therefore, the effective disorder imposed by impurities is below the critical value (see Appendix B), above which the first-order phase transition between isotropic (or paranematic) and nematic phase is replaced by a gradual evolution of ordering on varying temperature.Furthermore, one sees that, for {, } = {1, 0.1}, the hysteresis is almost totally suppressed indicating gradual evolution of nematic ordering on varying temperature [17].Therefore, for this choice of parameters, the effective disorder is supercritical.
Next, we consider impact of  on nematic ordering.We first demonstrate impact of  on range of ordering.In Figure 2, we plot () dependencies using parameter set {, , } = {0.1, 3, 0.5}, for (i)  = 0, (ii)  = 0.5, and (iii)  = 0.5 and both histories FQH and TQH.The parameter set {, , } was chosen in such a way that SRO is realized for (i)  = 0.In the case (ii) the external field is strong enough to replace SRO with QLRO.For still stronger fields, see (iii); LRO is realized.The range of order was obtained via finite size analysis, using (5).For example, in Figure 3(a), we show the finite size analysis for TQH.For (i)  = 0 (Figure 3(a)), we obtain  = 1.54±0.05fingerprinting SRO.For (ii)  = 0.5, (Figure 3(b)) simulations yield  = 0.34 ± 0.05, suggesting QLRO.On the other hand for (iii)  = 1.5, we get  = 0.00 ± 0.02, indicating LRO (Figure 3(c)).In Table 1 we summarize simulation results for  on varying  for extreme histories TQH and FQH.
Next, we consider impact of  on characteristic linear length  of the system.The following behavior is expected.For relatively weak external field strengths, the impurity imposed disorder dominantly influences nematic structure.According to Larkin-Imry-Ma theorem one expects a domain pattern characterized by   ∼  −2/(4−) , where  measures the disorder strength.Our mean-field estimate shown in Appendix B suggests  ∝ √ (see (B.6)).On the contrary, for relatively large external fields, one expects that the external field coherence length   ∝ 1/ characterizing the impact of  in bulk LC system would be observed.Using our simulations we estimate  by fitting () dependences, calculated from fixed point configurations, with (4), where  ∼ 1/.Characteristic dependencies () are shown in Figure 4.One indeed observes tendency  ∼   ∝ 1/ ( ∼   , roughly independent of ) for relatively strong (weak) external field .We refer to regimes where  ∼      , corresponding to SRO-QLRO transition, disorderexternal field dominance crossover, and QLRO-LRO transition, respectively.With increasing temperature these lines are pushed towards larger values of .Namely, due to increased thermal fluctuations, higher external fields are needed to reach regime exhibiting qualitatively stronger degree of ordering.
Finally, we analyzed external field controlled memory effects in systems.We performed the following -loops.(i) Initial configurations were obtained for a given parameter set {, , } via TQH for  = 0. (ii) Then we switched on the external field of strength , for which we calculated a fixed point configuration.(iii) Afterwards we switched the field off and calculated the corresponding fixed point nematic structure.We refer to respective configurations as the (i) reference, (ii) field-on, and (iii) field-off structures, respectively.In Figure 6, we plot the average degree of orientational order  obtained from  () on increasing  for different anchoring strengths deep in the nematic phase for  = 3.For this anchoring strength the reference structures exhibit SRO, yielding  = 0 within numerical error.One sees that in fieldon state the degree of ordering is substantially enhanced due to external field induced preference alignment.In field-off structures finite degree of ordering is preserved providing  < 0.2.We refer to the resulting ordering as remanent order   .
In Figure 6 one sees that   () displays two different regimes.
In the weak field regime   () monotonously increases with  until the saturation regime is entered, in which   () variations are relatively small on varying .

Conclusions
We study numerically memory effects using modified Random Anisotropy Nematic-type lattice model [16,17].In our model impurities are randomly spatially distributed enforcing to nearby nematic spins quenched randomly chosen orientation.In simulations we varied concentration  of impurities, coupling strength  between nematic spins and impurities, temperature, and external field strength .Furthermore, in most simulations we used two extreme histories of nematic spins, namely, field quenched history (FQH) and temperature quenched history (TQH).Using them, we gain information on memory capabilities of systems of interest.We calculated fixed point nematic configurations, revealing interplay between nematic interactions, local disordering, and global ordering fields.From fixed point configurations we calculated two-point orientational correlation function () and average local and global nematic order parameter.From calculated () profiles we determined characteristic length  of the system.Furthermore, using finite size analysis, we determine range of orientational ordering for a given set of parameters.
Key results of our study are as follows.We determined regimes exhibiting LRO, QLRO, or SRO orders as a function of temperature and external ordering field .The range of order was determined from values of  using finite size analysis.Note that the Lakrin-Imry-Ma [1,6] theorem predicts SRO in absence of .The theorem was derived for random-field type disorder.However, we used slightly different type of disorder, which is closer to random-bond type imposed frustrations.In simulations we obtained SRO only for strong enough disorder.Furthermore, our simulations confirm that range of order strongly depends on history of samples.In particular, we demonstrate that strong enough applied field  could change range of order.In the {, } plane we estimate critical line separating SRO-QLRO and QLRO-LRO for different histories of samples.Histories play significant role revealing glass-type properties of the systems.Furthermore, we determined impact of  on characteristic linear correlation length  within samples.We distinguish between disorder dominated and external field dominated regime.In the first regime impact of  is negligible.On the other hand in the second regime  is roughly proportional with 1/.We roughly determined crossover region separating these regimes for different temperatures.We also study impact of temporarily present field  on remanent average nematic ordering   .For this purpose we used parameter for which   = 0 (within numerical order) in the reference state obtained in absence of .We determined magnitude of   as a function of .Our simulations reveal that   monotonously increases with  until saturation region is reached.
Displayed memory capabilities could be of use for various external field driven memory devices.Modelled systems mimic to some extent LCs confined to various porous matrices or mixtures of nematic LCs with appropriate nanoparticles (e.g., aerosils).In our future studies we intend to explore in more detail analogies between our model and random magnets and spin glasses.

A. Brownian Dynamics-Type Numerical Simulation
We consider the following Hamiltonian: where  is given by ( 1) and   are Lagrange multipliers taking into account the normalization constraint | ⃗   | = 1.Within a sweep the orientations of all nematic spins are updated using the following equations: Here ⃗  (old)  and ⃗  (new)  represent orientations of nematic spins before and after the sweep, Δ measures the time step of each sweep,  is the representative LC rotational diffusion constant,   is the Boltzmann constant, is the so called residuum, and the quantity  ⇀  is defined in terms of unit vectors  ⇀ V 1 and  ⇀ V 2 as follows: The first term on the right side of (A.2) corresponds to a mechanical elastic torque tending to rotate the nematic spin towards a local energy minimum.The second term Δ ⃗  ()    mimics random thermal fluctuations.
In simulations we use the dimensionless time step Δ * and dimensionless temperature  * defined as follows: Note that when the residuum drops to zero (i.e., it is small enough with respect to equilibrium-like fluctuations), the system reaches a fixed point configuration reflecting either a metastable or globally stable nematic configuration.
To calculate a thermal fluctuation vector Δ ⃗  ()  , we first consider rotation of a nematic spin in its eigen frame, within which the local -axis points along the spin direction.Therefore, only rotations of the spin about perpendicular local and -axes are relevant and we set that both rotations are mutually independent.We express the corresponding rotational energy contribution as follows: Here Ω  and Ω  are the corresponding rotation angles and  1 is the appropriate dimensionless parameter which enables Advances in Condensed Matter Physics phenomenological tuning of thermal fluctuations.Assuming the canonical distribution, we express the infinitesimal probability for the range of angles within the phase-space element Ω  Ω  as follows: where  is the normalization constant.By expressing (A.6) in the two-dimensional polar coordinate system it follows Here Ω is the magnitude of the vector (Ω  , Ω  ) and the angle Ψ defines its direction with respect to the local -axis.Therefore, we mimic thermal fluctuations by the rotational kinetic free energy with two independent components of the angular velocity   and   .Within a sweep, each spin is rotated by the angle Ω about the axis that is perpendicular to the local -axis and makes the angle Ψ with the local -axis.
The distribution of Ψ is uniform in the interval (0, 2) and the probability distribution function for Ω is given by After the thermal rotation of the spin in its own frame is calculated the rotated spin is transformed to the laboratory system using the following transformation: (A.9)

B. Phenomenological Estimates
In the appendix we present estimates on phase and structural LC behavior in presence of impurity spins, where we assume that the system exhibits a domain-type pattern characterized by a single characteristic length.Randomly distributed impurities are assumed to produce random field-type disturbances.
Impurities are homogeneously distributed within the sample of volume  with the volume concentration  =  im V im /.Here V im determines the volume of an impurity and  im is the total impurity number inside .We estimate the nematic ordering at mesoscopic level by the uniaxial tensor order parameter  = ( ⃗  ⊗ ⃗  − /3).The quantity ⃗  is the nematic director field and  is the uniaxial orientational order parameter.The corresponding total free energy is given by The first integral is carried over the LC volume, and the second one over LC-impurity interfaces.
The condensation (  ), elastic (  ), and interfacial (  ) terms are expressed as follows:   | = 1.We assume that impurities are essentially homogeneously distributed with concentration  and that orientational probability distribution of ⃗   is spatially isotropic.We adopt the Larkin-Imry-Ma picture in which the disorder breaks the system into domains of characteristic size .Our aim is to express it in terms of system parameters.For this purpose we considered the mean free energy of a single domain of volume   ∼  3 , within which there are  ()  im impurities possessing surface area  im .Within each domain the nematic director is regarded to be on average aligned along a single symmetry breaking direction.The corresponding free energy penalty within a domain volume is estimated by Here ⟨ 2 ⟩ = (1/2)⟨3( ⃗  ⋅ ⃗   ) 2 − 1⟩ is the mean value of the second Legendre polynomial of ( ⃗ ⋅ ⃗   ) within the domain.The penultimate term evaluates approximately the effect of the gradient term.The final term in (B.3) comes from evaluating   over the domain.The final term can be expressed in terms of  as follows.We assume the isotropic orientational distribution of orientations ⃗   .In an infinitely large domain one expects ⟨ 2 ⟩ = 0.In a finite volume the cancellation is only partial.According to the central limit theorem one expects ⟨ 2 ⟩ ∼ ( ()  im ) 1/2 . The number of impurities and the characteristic impurity separation  im are related by  ()  im ∼ (/ im ) 3 , where  im ∼ (V im /) 1/3 .An expression for  is obtained by balancing the  -dependent elastic and interface contributions in Δ  .It follows where  = (√/3)( im /√V im )( (0)2  / (0)   3/2 ) and  eff =  +  (0)2  /2 2 play roles of dimensionless field and effective dimensionless temperature, respectively.In the limit of strong enough anchoring, where  ∼  im ∼ (V im /) 1/3 , and for spherical impurities of radius , one obtains It is now possible to link this model to the theory of nematics in a constant external field.According to it the 1st order phase transition exists at critical temperature corresponding to  eff = 1 +  as long as  <   ≡ 0.5.For  >   one obtains gradual evolution of ordering on varying temperature.

Figure 6 :
Figure 6: Average degree of ordering for switched-on and afterwards switched-off external ordering field  for different values of .Full symbols:  is switched on.Empty symbols:  is switched off.TQH,  = 3,  = 05, and  = 70.