format: text

text:
# Universality and scaling study of the critical behavior of the two-dimensional Blume-Capel model in short-time dynamics

## 1 INTRODUCTION

Numerical simulation in the short-time regime has become an important tool to study phase transitions and critical phenomena. The reason is that universality and scaling behavior are already present in the dynamic systems since the early stages of their evolution . Moreover, this kind of approach reveals the existence of a new and unsuspected critical exponent. As shown by Janssen et al. on basis of renormalization group theory, if we tune the parameters at their critical values but with initial configurations characterized by nonequilibrium states, the time evolution of quantities like magnetization exhibit a polynomial behavior governed by an exponent $`\theta `$, which is independent of the known set of static exponents and of the dynamical critical exponent $`z`$. This new exponent characterizes the so called “critical initial slip”, the anomalous behavior of the magnetization when the system is quenched to the critical temperature $`T_c`$. Working with systems without conserved quantities, model A in the terminology of Halperin et al. , Janssen et al. found a scaling form for the moments of the magnetization, which sets soon after a microscopic time scale $`t_{mic}`$. Those relations have been confirmed in several numerical experiments . For the $`k\mathrm{th}`$ moment of the magnetization, this scaling form reads

$$M^{(k)}(t,\tau ,L,m_0)=b^{k\beta /\nu }M^{(k)}(b^zt,b^{1/\nu }\tau ,b^1L,b^{x_0}m_0).$$
(1)
Here $`b`$ is an arbitrary spatial scaling factor, $`t`$ is the time evolution and $`\tau `$ is the reduced temperature, $`\tau =(TT_c)/T_c`$. As usual, the exponents $`\beta `$ and $`\nu `$ are the well-known static exponents, whereas $`z`$ is the dynamic one. Equation (1) depends on the initial magnetization $`m_0`$ and gives origin to the new exponent $`x_0`$, the scaling dimension of the initial magnetization, related to $`\theta `$ by $`\theta =(x_0\beta /\nu )/z`$.

For a large lattice size $`L`$ and small initial magnetization $`m_0`$, the system in its early stage presents small spatial and temporal correlation lengths, which may eliminate usual finite size problems. In this limit, if we choose the scaling factor $`b=t^{1/z}`$  at the critical temperature ($`\tau =0`$), we obtain

$$M(t,m_0)m_0t^\theta $$
(2)
from the scaling relation (1). The exponent $`\theta `$ has been calculated for the two-dimensional (2d)  and three-dimensional (3d)  Ising models, 2d 3-state Potts model , Ising model with next-nearest neighbor interactions  and Ising model with a line of defects . In addition, this short-time universal behavior was found in irreversible models with synchronous  and continuous time dynamics . In all of those cases, a positive value for $`\theta `$ has been found, that indicates a surprising initial increasing of the magnetization in the short-time regime $`t_0m_0^{z/x_0}`$. This effect can be related to a “mean field” behavior since the system presents small correlation length in the beginning of the time evolution. Thus, when the system is quenched to the critical temperature $`T_c`$, it feels as being in an ordered state since $`T_c<T_c^{MF}`$ .

On the other hand, as shown by Janssen and Oerding , the behavior of a thermodynamic system is more complex at a tricritical point (TP) and the corresponding exponent $`\theta `$ may attain negative values.

At a tricritical point the magnetization shows a crossover from the logarithmic behavior $`M(t)m_0(\mathrm{ln}(t/t_0))^a`$ at short times $`t<<m_0^4`$ to $`t^{1/4}`$ power law with logarithmic corrections, $`M(t)(t/\mathrm{ln}(t/t_0))^{1/4}`$ in 3 dimensions. The above behavior can be stated in the generalized form

$$M(t)=m_0(\mathrm{ln}(t/t_0))^aF_M\left(\left(\frac{t}{\mathrm{ln}(t/t_0)}\right)^{1/4}(\mathrm{ln}(t/t_0))^am_0\right),$$
(3)
where $`F_M(x)1`$ or $`F_M(x)1/x`$, respectively for vanishing and large arguments. Below 3 dimensions it reduces to the scaling form

$$M(t)m_0t^\theta .$$
(4)
Here $`\theta `$ is the exponent related to the tricritical point of the relaxation process at early times and it is expected to assume negative values.

In this paper, we perform short-time Monte Carlo (MC) simulations to explore the critical dynamics of the 2d Blume-Capel model. We evaluate the dynamic exponents $`\theta `$ and $`z`$, besides the static exponents $`\nu `$ and $`\beta `$ at the tricritical point. To the best of our knowledge this is the first time it is done numerically. We also estimate the dynamic exponents along the critical line. We observe a clear trend toward the values of $`z`$ and $`\theta `$ for the corresponding 2d Ising values when the crystal field $`D`$ becomes large and negative, indicating dynamic universality in the limit $`D\mathrm{}`$.

The paper is organized as follows. In the next section we present the model and its phase diagram. Sec. III contains the main scaling relations and describes our short-time MC simulations. Results are presented for critical points on the second-order line. In Sec. IV, we explore the short-time dynamics to study the tricritical behavior. Sec. V contains a brief outlook and concluding remarks.

## 2 The model

The Blume-Capel  (BC) model is a spin-1 model which has been used to describe the behavior of $`{}_{}{}^{3}\mathrm{He}{}_{}{}^{4}\mathrm{He}`$ mixtures along the $`\lambda `$ line and near the critical mixing point. Apart from its practical interest, the BC model has intrinsic interest since it is the simplest generalization of the Ising model ($`s=1/2`$) exhibiting a rich phase diagram with first and second-order transition lines besides a tricritical point. In real systems tricritical points appear in $`{}_{}{}^{3}\mathrm{He}{}_{}{}^{4}\mathrm{He}`$ mixtures, such that when a small fraction of $`{}_{}{}^{3}\mathrm{He}`$ is added to $`{}_{}{}^{4}\mathrm{He}`$, a critical line terminates at a concentration of $`{}_{}{}^{3}\mathrm{He}`$ approximately at $`0.67`$. The BC model or its well known generalization, the Blume-Emery-Grifftihs model  was studied by mean-field approximation, real space and renormalization group schemes , Monte Carlo renormalization group approach  and finite-size scaling combined with conformal invariance . The Hamiltonian of the two-dimensional model is

$$H=J\underset{<i,j>}{}S_iS_j+D\underset{i=1}{}S_i^2,$$
(5)
where $`<i,j>`$ indicates nearest neighbors on $`L^2`$ lattices and $`S_i=\{1,0,1\}`$. The parameter $`J`$ is the exchange coupling constant and $`D`$ is the crystal field. We show its phase diagram in Fig. 1. Table I lists points on the second order critical line and the tricritical point where we have performed our simulations. Points in Table I were obtained from  and from a private communication of the authors in . That table also contains our results for the corresponding critical and tricritical exponents.

We remark that along the critical line, this model presents a critical behavior similar to that of the Ising model. However, exactly at the tricritical point the exponents change abruptly. They are given by the dimensions of the irreducible representations of the Virasoro algebra  with central charge (conformal anomaly number) $`c=7/10`$ . In  finite-size scaling combined with conformal invariance permited to observe a smooth change between Ising-like and tricritical behavior. In finite systems, Ising-like behavior is reached only when $`D\mathrm{}`$. In that limit $`\eta /2(21/\nu )0.125`$ that is the exact value for the Ising model. In our short-time simulations the same kind of crossover behavior is observed for the dynamic exponents $`z`$ and $`\theta `$ when we move along the critical line.

## 3 Non-equilibrium short-time dynamics at a critical point

In short-time MC simulations critical slowing down can be neglected. It happens because spatial and time correlation lengths are small in the early stages of evolution. On the other hand, we need to deal with several samples taken over independent initial configurations since the systems which are being simulated are far from equilibrium. In fact this approach requires calculation of average (over samples) magnetization and of its moments $`M^{(k)}(t)`$,

$$M^{(k)}(t)=\frac{1}{N_sL^{kd}}\underset{j=1}{\overset{N_s}{}}\left(\underset{i=1}{\overset{L^d}{}}\sigma _{ij}(t)\right)^k,$$
(6)
where $`\sigma _{ij}(t)`$ denotes the value of spin $`i`$ of $`j`$th sample at the $`t`$th MC sweep. Here $`N_s`$ denotes the number of samples and $`L^d`$ is the volume of the system. This kind of simulation is performed $`N_B`$ times to obtain our final estimates in function of $`t`$. In this paper, the dynamic evolution of the spins $`\{\sigma _i\}`$ is local and updated by the heat-bath algorithm.

### A The critical initial slip

The evolution of the $`k`$th moment of magnetization in the initial stage of the dynamic relaxation can be obtained from Eq. (1) for large lattice sizes $`L`$ at $`\tau =0`$ with $`b=t^{1/z}`$. This yields

$$M^{(k)}(t,m_0)=t^{k\beta /\nu z}M^{(k)}(1,t^{x_0/z}m_0).$$
(7)
By expanding the corresponding first moment equation for small $`m_0`$, we obtain Eq. (2) under the condition that $`t^{x_0/z}m_0`$ is sufficiently small, which sets a time scale $`t_0m_0^{z/x_0}`$  where that phenomena can be observed.

In Fig. 2a we present our results for the exponent $`\theta `$ at the critical point $`k_BT_c/J=1.695`$ and $`D_c/J=0`$, for lattice size $`L=80`$ and 5 different initial magnetizations $`m_0`$. Our estimates for each $`\theta =\theta (m_0)`$ were obtained from $`N_B=5`$ independent bins with $`N_S=10000`$, for $`t`$ up to 100 sweeps. Figure 2b illustrates the numerical evaluation of $`\theta `$ for $`m_0=0.02`$ from a log-log plot of the magnetization versus time. The linear fitting in Fig. 2a gives $`\theta =0.193(2)`$ with goodness of fit  $`Q=0.72`$.

Another method has been recently proposed by Tomé and de Oliveira  to evaluate $`\theta `$. It avoids the sharp preparation of samples with defined and nonzero magnetization and the delicate numerical extrapolation $`m_00`$. The method is based on the time correlation function of the total magnetization,

$$C(t)=\frac{1}{L^d}\underset{i=1}{\overset{L^d}{}}\underset{j=1}{\overset{L^d}{}}\sigma _i(t)\sigma _j(0).$$
(8)
Starting from random initial configurations the above correlation behaves as $`C(t)t^\theta `$, which permits us to obtain the exponent $`\theta `$ from a log-log plot of $`C(t)`$ versus $`t`$. We obtained $`\theta =0.194(3)`$ for $`k_BT_c/J=1.695`$ and $`D_c/J=0`$ choosing the time interval $`[20150]`$ due to the highest value of $`Q`$ ($`Q=0.99`$). This value is in complete agreement with our above estimate of the exponent $`\theta `$ and it is consistent within error bars with previous results for the 2d Ising model. In table I we also present results for $`\theta `$ at several points of the critical line.

### B Dynamic critical exponent $`z`$

The observables in short-time analysis are described by different scaling relations according to the initial magnetizations. In particular, the second moment $`M^{(2)}(t,L)`$ in Eq. (6),

$$M^{(2)}=\frac{1}{L^{2d}}\underset{i=1}{\overset{L^d}{}}\sigma _i^2+\frac{1}{L^{2d}}\underset{ij}{\overset{L^d}{}}\sigma _i\sigma _j,$$
(9)
with $`m_0=0`$ behaves as $`L^d`$ since in the short-time evolution the spatial correlation length is very small when compared with the lattice size $`L`$. Thus, we arrive at

$$M^{(2)}(t,L)=t^{2\beta /\nu z}M^{(2)}(1,t^{1/z}L)t^{(d2\beta /\nu )/z}.$$
(10)
This equation can be used to determine relations involving static critical exponents and the dynamic exponent $`z`$ . However, a way to evaluate independently the exponent $`z`$ is through out the time-dependent fourth-order Binder cumulant at the critical temperature ($`\tau =0`$),

$$U_4(t,L,m_0)=1\frac{M^{(4)}(t,L,m_0)}{3\left(M^{(2)}(t,L,m_0)\right)^2},$$
(11)
which obeys the equation

$$U_4(t,L,m_0)=U_4(b^zt,b^1L,b^{x_0}m_0).$$
(12)
If we set $`m_0=0`$, we eliminate the dependence on the exponent $`x_0`$ and the exponent $`z`$ can be evaluated through scaling collapses of the generalized cumulant for different lattice sizes . To match the Binder cumulants $`U_4(t_1,L_1)`$ and $`U_4(t_2,L_2)`$ obtained from two time series for lattice sizes $`L_1`$ and $`L_2`$, respectively, with $`b=L_2/L_1`$ ($`L_2>L_1`$), we interpolate the series $`U_4(t,L_1)`$ to obtain $`\stackrel{~}{U}_4(b^zt,L_1)`$. Next, we define the function

$$\chi ^2(z)=\frac{1}{t_ft_i}\underset{t=t_i}{\overset{t_f}{}}\left[\stackrel{~}{U}_4(b^zt,L_1)U_4(t,L_2)\right]^2,$$
(13)
where the best estimate for $`z`$ corresponds to the one which minimizes $`\chi ^2(z)`$.

In Fig. 3 we show the scaling collapses of the Binder cumulants for different lattice sizes. We have collapsed the following pairs of lattices $`(L_1,L_2)=(10,20),(20,40)`$ and $`(40,80)`$. From the largest pair of lattices we obtained $`z=2.159(6)`$ in the time interval $`[501000]`$. Our final error estimate is based on 25 different collapses obtained from $`N_B=5`$ independent bins for each lattice size.

Another universal behavior of the dynamic relaxation process also described by Eq. (1) can be obtained with the initial condition $`m_0=1`$ . This condition is related to another fixed point in the context of renormalization group approach. Thus, starting from an initial ordered state one obtains a power law decay of the magnetization at the critical temperature,

$$M(t)t^{\beta /\nu z},$$
(14)
when we choose $`b^zt=1`$ in the limit of $`L\mathrm{}`$. Taking into account this relation, another method has been proposed  to estimate the dynamic exponent $`z`$. This approach uses the second cumulant

$$U_2(t,L)=\frac{M^{(2)}(t,L)}{(M(t,L))^2}1,$$
(15)
which should take the simple form $`U_2(t,L\mathrm{})t^{d/z}`$. The advantage of this procedure is that curves for different lattices lay on the same straight line in a log-log plot without any re-scaling in time. However, this technique has not been successful in at least two well known models: the two-dimensional $`q=3`$ Potts model  and the Ising model with three spin interactions in just one direction . The reason for the above disagreement may be related to the behavior of the second term of r.h.s in Eq. (9) when $`m_0=1`$. We have proposed  that this behavior could indeed be obtained working with the ratio $`F_2=M^{(2)}/M^2`$ with different initial conditions for each moment. As we know the behavior of the second moment of the magnetization when samples are initially disordered $`(m_0=0)`$ and also the time dependence of the magnetization when samples are initially ordered $`(m_0=1)`$, we easily obtain

$$F_2(t)=\frac{M^{(2)}(t,L)|_{m_0=0}}{(M(t,L))^2|_{m_0=1}}t^{d/z}.$$
(16)
A log-log plot with error bars for the critical point $`k_BT_c/J=1.695`$ and $`D_c/J=0`$ is presented in Fig. 4 for $`L=160`$. We obtained $`z=2.1057(7)`$ with $`Q=0.99`$ in the range $`[30200]`$, which does not agree with the value obtained from Eq. (12). However, as we move away from the tricritical point, the values of $`z`$ obtained (Table 1) with Eq. (16) show a clear trend toward the expected value of the dynamic exponent $`z`$ $`(z=2.1567(7))`$ of the 2d Ising model . On the other hand, the values obtained from the cumulant in Eq. (12) remain essentially the same along the entire critical line. It seems the cumulant $`U_4`$ is less sensitive to crossover effects than the mixed method.

### C Static exponents and universality class

The exponent $`1/\nu z`$ can be obtained by differentiating $`\mathrm{ln}M(t,\tau ,m_0)`$ with respect to the temperature at $`T_c`$,

$$\frac{\mathrm{ln}M(t,\tau ,L)}{\tau }|_{\tau =0}t^{1/\nu z},$$
(17)
if we consider the scaling relation for the magnetization when the initial state of samples is ordered ($`m_0=1`$) .

Our results for $`\nu z`$ were obtained through finite diferences at $`T_c\pm \delta `$ with $`\delta =0.001`$. They rely on $`N_B=5`$ independent bins with $`N_s=5000`$ samples each for $`L=160`$. The time interval $`[80200]`$ corresponds to the range where the goodness of fit parameter attains its highest value $`(Q=0.99)`$.

In Fig. 5 we show the log-log behavior of the derivative $`_\tau \mathrm{ln}M(t)`$ at $`k_BT_c/J=1.695`$ and $`D_c/J=0`$. In Table I we present our final estimates of $`1/\nu `$. The 6th column is obtained with the estimates of $`z`$ from Eq. (12) (data in 4th column), while the 7th column corresponds to estimates for $`\nu `$ with values of $`z`$ from Eq. (16) (data in 5th column).

Since we have already collected estimates for $`\nu z`$, it is straightforward to obtain estimates for $`\beta `$ following Eq. (14). Our estimates of $`\beta `$ are presented in the last column. Our values in Table I can be compared with theoretical predictions for an Ising like critical point ($`1/\nu =1,\beta =1/8`$).

## 4 Results from short-time dynamics at the tricritical point

From the results presented in Refs.  and  we can describe the time dependence of the magnetization $`(k=1)`$ for the 2d Blume-Capel model as

$$M(t)\{\begin{array}{ccccc}m_0t^\theta ,& 0<t<t_0,& \mathrm{where}\theta 0& & \\ t^{\beta /\nu z},& t_0<t<t_I& & & \\ & & & & \\ m_0t^\theta ,& 0<t<t_0,& \mathrm{where}\theta 0& & \\ t^{\beta /\nu z},& t_0<t<t_I& & & \end{array}\begin{array}{c}\text{critical}\\ \text{point}\\ \\ \text{tricritical}\\ \text{point}\end{array}$$
(18)
for an initial small magnetization $`m_0`$. Here $`t_I`$ stands for the time before the system has reached the thermalization.

We also included in Table I our estimates for $`\theta ,z,1/\nu `$ and $`\beta `$ at the tricritical point $`k_BT_t/J=0.610`$ and $`D_t/J=1.9655`$ working with lattice size $`L=80`$.

In Fig. 6a we show the values of $`\theta `$ for 5 different initial magnetizations $`m_0`$ at the tricritical point and in Fig. 6b the time evolution of the magnetization for $`m_0=0.08`$ in a log-log plot. Our estimates for each $`\theta (m_0)`$ were obtained from $`N_B=20`$ independent bins with $`N_S=10000`$ samples, for $`t`$ up to 80 sweeps. The linear extrapolation in Fig. 6a gives $`\theta =0.53(2)`$ with goodness of fit $`Q=0.75`$. The corresponding study with the time correlation function in Eq. (8) also gives $`\theta =0.53(2)`$ with $`Q=0.99`$ in the time interval $`[2080]`$. The generalization of the dynamic scaling relation for the $`k`$th moment of the magnetization at a tricritical point can be written as

$$M^{(k)}(t,\tau ,g,L,m_0)=b^{k\beta _t/\nu }M^{(k)}(b^zt,b^{1/\nu }\tau ,b^{\varphi _t/\nu }g,b^1L,b^{x_0}m_0).$$
(19)
It differs from the critical case by the scaling field $`g`$ that measures the shunting line of critical point along of tangent transition line (at the tricritical point). The quantity $`\varphi _t`$ is the known as the crossover exponent. At the criticality we have $`\varphi _t=g=0`$.

We show in Fig.7 scaling collapses of the cumulant $`U_4(t,L)`$ at the tricritical point, whose behavior is quite different of the corresponding one at a critical point (Fig. 3). Our first estimate based on Eq. (12) leads to $`z=2.21(2)`$. This value was obtained from the pair of largest lattice sizes $`L=40`$ and $`80`$. The value remains the same when we consider the time interval $`[10,1000]`$ or $`[200,1000]`$. The second estimate for the dynamic exponent, based on Eq. (16) gives $`z=2.215(2)`$. It was obtained from a larger lattice ($`L=160`$) in the time interval $`[30,200]`$, with $`Q=0.71`$. We do not show the log-log plot of $`F_2(t)`$ in this case because it is quite similar to Fig. 4. We had to restrict the time interval, when compared with the $`U_4`$ calculation, in order to obtain acceptable values for $`Q`$. Here, in contrast to the estimates of $`z`$ (4th and 5th column) at the critical points listed in Table I, both methods lead to the same numerical result.

The evaluation of static exponents $`\nu `$ and $`\beta `$ at the tricritical point follows the same procedure as applied to the critical line. Estimates exhibited in Table I are in good agreement with results provided by conformal invariance $`\nu =9/5`$ and $`\beta =1/24`$.

Next, in order to study the influence of the local dynamics over the values of the exponents we recall the simulations with Glauber dynamics performed by de Alcantara Bonfim. Our estimates of $`\beta /\nu z`$ obtained from the decay of the magnetization (starting from initialy ordered state) is $`0.0381(1)`$ with the heat-bath dynamics while $`\beta /\nu z=0.0377(7)`$ when updated following Glauber rules.

## 5 CONCLUSIONS

We have performed short-time Monte Carlo simulations to evaluate dynamic and static exponents at critical and tricritical points of the spin-1 Blume-Capel model. According to analytical predictions by Janssen and Oerding, a negative value for the new exponent $`\theta `$ was obtained for the tricritical point. The dynamic exponent $`z`$ was estimated by colapsing the fourth-order Binder cumulant and also by the function $`F_2(t)`$ which explores time evolution of the magnetization starting from different initial conditions. The value of $`z`$ was used to obtain the static exponents $`\beta `$ and $`\nu `$ which resulted in good agreement with the exact results provided by the tricritical Ising model . The dynamic exponents also were calculated along the critical line. Estimates for $`\theta `$ are in good agreement with known result for the 2d Ising model. However, near the tricritical point we found strong crossover effects on dynamical exponent $`z`$ when using the recently proposed  technique based on mixed initial conditions. On the other hand, the fourth-order Binder cumulant is almost free from crossover effects.

Acknowledgments

R. da Silva and J. R. Drugowich de Felício gratefully acknowledge support by FAPESP (Brazil) and N.A. Alves by CNPq (Brazil). The authors deeply thank also to J. A. Plascak and J. C. Xavier for their communication about the values of critical points for $`D<0`$ and DFMA (IFUSP) for the computer facilities extended to us.


text:
# References

RH-03-2002

A Note on Noncommutative Scalar Multisolitons

Bergfinnur Durhuusa111email: durhuus@math.ku.dk and Thordur Jonssonb222e-mail: thjons@raunvis.hi.is

aMatematisk Institut, Universitetsparken 5

2100 Copenhagen Ø, Denmark

bUniversity of Iceland, Dunhaga 3,

107 Reykjavik, Iceland

Abstract. We prove that there do not exist multisoliton solutions of noncommutative scalar field theory in the Moyal plane which interpolate smoothly between $`n`$ overlapping solitons and $`n`$ solitons with an infinite separation.

Solitons in scalar field theories were first studied in  in the limit of an infinite noncommutativity parameter $`\theta `$ where the kinetic term in the action can be neglected. The existence of rotationally invariant scalar noncommutative solitons at finite but sufficiently large $`\theta `$ was proved in  and stability of these solitons was established in . In  the moduli space of the infinite $`\theta `$ solitons is studied and it is argued that solitons attract at finite $`\theta `$ so that all multi-soliton solutions at finite $`\theta `$ are rotationally invariant around some point. In this paper we place this conjecture on a firm foundation. More precisely, we show that there does not exist a family of solitonic solutions interpolating smoothly between two overlapping solitons and two infinitely separated ones. The proof is based on using energy estimates which indeed show that the energy of two infinitely separated solitons is greater than that of two overlapping ones. Various aspects of the theory of solitons in noncommutative scalar field theories are discussed in .

Solitons in a noncommutative two-dimensional scalar field theory with a potential $`V`$ are finite energy solutions to the variational equations of the energy functional

$$S(\phi )=\mathrm{Tr}([a,\phi ][\phi ,a^{}])+\theta V(\phi )),$$
(1)
where $`a^{}`$ and $`a`$ are the usual raising and lowering operators of the simple harmonic oscillator and $`\phi `$ is a self-adjoint operator on $`L^2(𝐑^2)`$. We assume as in  that the potential $`V`$ is at least twice continuously differentiable with a second order zero at $`x=0`$, $`V(x)>0`$ if $`x0`$, and that $`V`$ has only one local minimum in addition to the one at $`x=0`$ located at $`s>0`$. For technical reasons we assume that the potential is analytic on a neighbourhood of the interval $`[0,s]`$ but we do not believe that this condition is necessary. We follow the notation of .

Let $`\varphi _1`$ and and $`\varphi _2`$ denote the stable one and two soliton solutions which were constructed for sufficiently large $`\theta `$ in . These solitons have the properties

$$\varphi _1s|00|,\varphi _2s(|00|+|11|),$$
(2)
as $`\theta \mathrm{}`$, i.e., they converge to the infinite $`\theta `$ solitons of . Let $`U_z`$ be the unitary operator on $`L^2(𝐑^2)`$ which implements the translation of the origin in the plane to the point $`(x,y)`$ and $`z=x+iy`$. Then $`|z=U_z|0`$ is a usual harmonic oscillator coherent state and the operator $`s|zz|`$ is interpreted as an infinte $`\theta `$ soliton located at $`z`$. Similarly, one can argue  that

$$\varphi _{2,z}^{(\mathrm{})}=s(|z_+z_+|+|z_{}z_{}|),$$
(3)
where

$$|z_\pm =\frac{|z\pm |0}{\sqrt{2(1\pm e^{|z|^2/2})}},$$
(4)
is a two soliton solution at infinte $`\theta `$ with the solitons located at $`z`$ and the origin. Clearly

$$(\varphi _{2,z}^{(\mathrm{})}s(|zz|+s|00|)0\mathrm{as}|z|\mathrm{}$$
(5)
and

$$\varphi _{2,z}^{(\mathrm{})}s(|00|+|11|)\mathrm{as}|z|0.$$
(6)
It is natural to interpret the operator

$$\varphi _{1,z}=U_z\varphi _1U_z^{}$$
(7)
as the finite $`\theta `$ one soliton solution located at $`z`$. The question is whether there is a finite $`\theta `$ analogue of $`\varphi _{2,z}^{(\mathrm{})}`$. In  it is argued that this is not the case. The argument is based on computing the energy of the infinite $`\theta `$ solution using the finite $`\theta `$ energy functional and observing that the energy depends on the separation between the solitons. Here we provide a proof.

Let us assume that there exists a family of solitons, $`\varphi _{2,x}`$, $`x0`$, depending smoothly on $`x`$, such that $`\varphi _{2,0}=\varphi _2`$ and

$$\underset{x\mathrm{}}{lim}(\varphi _{2,x}\varphi _{1,x}\varphi _1)=0,$$
(8)
where the convergence above is in the norm $`_{2,2}`$ defined by

$$\varphi _{2,2}^2=\mathrm{Tr}([a,\varphi ][\varphi ,a^{}]+\varphi ^2).$$
(9)
Then the energy $`E(x)=S(\varphi _{2,x})`$ is a differentiable function of $`x`$ and since $`\varphi _{2,x}`$ is a critical point of $`S`$ we have $`E^{}(x)=0`$ so $`E(x)`$ is a constant. Hence,

$$S(\varphi _2)=E(0)=\underset{x\mathrm{}}{lim}E(x)=2S(\varphi _1).$$
(10)
In order to prove the last equality in (10) we write

$$S(\varphi )=K(\varphi )+\theta \mathrm{Tr}V(\varphi ),$$
(11)
where

$$K(\varphi )=\mathrm{Tr}([a,\varphi ][\varphi ,a^{}].$$
(12)
We first show that

$$K(\varphi _{2,x})2K(\varphi _1)$$
(13)
as $`x\mathrm{}`$. By the triangle inequality

$`|K(\varphi _{2,x})^{1/2}K(\varphi _{1,x}+\varphi _1)^{1/2}|`$ $``$ $`K(\varphi _{2,x}\varphi _{1,x}\varphi _1)^{1/2}`$ (14)
$``$ $`\varphi _{2,x}\varphi _{1,x}\varphi _1_{2,2}.`$
In order to establish (13) it therefore suffices to show that $`K(\varphi _{1,x}+\varphi _1)2K(\varphi _1)`$ as $`x\mathrm{}`$. We have

$`K(\varphi _{1,x}+\varphi _1)`$ $`=`$ $`K(\varphi _{1,x})+K(\varphi _1)+2\mathrm{R}\mathrm{e}\mathrm{Tr}[a,\varphi _{1,x}][\varphi _1,a]`$ (15)
$`=`$ $`2K(\varphi _1)+2\mathrm{R}\mathrm{e}\mathrm{Tr}[a,\varphi _{1,x}][\varphi _1,a].`$
If $`f`$ is the function on $`𝐑^2`$ which corresponds to $`\varphi _1`$ under Weyl quantization, then

$$\mathrm{Tr}[a,\varphi _{1,x}][\varphi _1,a]=\frac{1}{4\pi }_{𝐑^2}_{\overline{z}}f(zx)_z\overline{f}(z)dzd\overline{z}$$
(16)
and the last integral tends to zero as $`x\mathrm{}`$ by the Riemann-Lebesgue Lemma (this can also be seen directly since $`f`$ is a gaussian).

We next prove the convergence

$$V(\varphi _{2,x})2V(\varphi _1)$$
(17)
as $`x\mathrm{}`$ for $`V`$ a polynomial. The proof generalizes without difficulty to the case where $`V`$ is analytic on a neighbourhood of $`[0,s]`$ by means of the analytic functional calculus. Let $`_2`$ denote the usual Hilbert-Schmidt norm. For $`n2`$ we have

$`|\mathrm{Tr}\varphi _{2,x}^n\mathrm{Tr}(\varphi _{1,x}+\varphi _1)^n|`$ $`=`$ $`\mathrm{Tr}\left({\displaystyle \underset{k=0}{\overset{n1}{}}}\varphi _{2,x}^k(\varphi _{2,x}\varphi _{1,x}\varphi _1)(\varphi _{1,x}+\varphi _1)^{n1k}\right)`$ (18)
$``$ $`{\displaystyle \underset{k=0}{\overset{n1}{}}}\varphi _{2,x}\varphi _{1,x}\varphi _1_2\varphi _{2,x}_2^k\varphi _{1,x}+\varphi _1_2^{n1k}`$
$``$ $`c\varphi _{2,x}\varphi _{1,x}\varphi _1_2,`$
for an $`x`$-indpendent constant $`c`$ since

$$\varphi _{x,1}+\varphi _1_2\varphi _{1,x}_2+\varphi _1_2=2\varphi _1_2.$$
(19)
On the other hand,

$$\mathrm{Tr}(\varphi _{1,x}+\varphi _1)^n2\mathrm{T}\mathrm{r}\varphi _1^n$$
(20)
as $`x\mathrm{}`$ since all cross terms between $`\varphi _{1,x}`$ and $`\varphi _1`$ vanish by the Riemann-Lebesgue Lemma as before. We have therefore established the convergence (17) for any polynomial $`V`$ with a second order zero at the origin and this completes our proof of (10). Now we return to our main line of argument.

Differentiating the equation $`S(\varphi _2)=2S(\varphi _1)`$ with respect to $`\theta `$ we obtain

$$\mathrm{Tr}V(\varphi _2)=2\mathrm{Tr}V(\varphi _1).$$
(21)
We claim that the above equality is violated for all sufficiently large values of $`\theta `$ and therefore the family $`\varphi _{2,x}`$ does not exist for large values of $`\theta `$.

Let $`\lambda _0^{(1)},\lambda _1^{(1)},\lambda _2^{(1)},\mathrm{}`$ and $`\lambda _0^{(2)},\lambda _1^{(2)},\lambda _2^{(2)},\mathrm{}`$ denote the eigenvalues of $`\varphi _1`$ and $`\varphi _2`$ in the harmonic oscillator basis. Then the eigenvalues form decreasing sequences of positive numbers and in  it is proven that

$`\lambda _0^{(1)}`$ $`=`$ $`s{\displaystyle \frac{2s}{\theta V^{\prime \prime }(s)}}+O(\theta ^2)`$ (22)
$`\lambda _1^{(1)}`$ $`=`$ $`{\displaystyle \frac{2s}{\theta V^{\prime \prime }(0)}}+O(\theta ^2)`$ (23)
$`\lambda _0^{(2)}`$ $`=`$ $`sO(\theta ^2)`$ (24)
$`\lambda _1^{(2)}`$ $`=`$ $`s{\displaystyle \frac{4s}{\theta V^{\prime \prime }(s)}}+O(\theta ^2)`$ (25)
$`\lambda _2^{(2)}`$ $`=`$ $`{\displaystyle \frac{4s}{\theta V^{\prime \prime }(0)}}+O(\theta ^2)`$ (26)
and all the other eigenvalues are $`O(\theta ^2)`$. It follows that

$`\mathrm{Tr}V(\varphi _1)`$ $`=`$ $`V(\lambda _0^{(1)})+V(\lambda _1^{(1)})+{\displaystyle \underset{n=2}{\overset{\mathrm{}}{}}}V(\lambda _n^{(1)})`$ (27)
$`=V(s)+\left({\displaystyle \frac{2s^2}{V^{\prime \prime }(s)}}+{\displaystyle \frac{2s^2}{V^{\prime \prime }(0)}}\right){\displaystyle \frac{1}{\theta ^2}}+{\displaystyle \underset{n=2}{\overset{\mathrm{}}{}}}V(\lambda _n^{(1)})+O(\theta ^3).`$ (28)
Similarly,

$`\mathrm{Tr}V(\varphi _2)`$ $`=`$ $`V(\lambda _0^{(2)})+V(\lambda _1^{(2)})+V(\lambda _2^{(2)})+{\displaystyle \underset{n=3}{\overset{\mathrm{}}{}}}V(\lambda _n^{(2)})`$ (29)
$`=2V(s)+\left({\displaystyle \frac{8s^2}{V^{\prime \prime }(s)}}+{\displaystyle \frac{8s^2}{V^{\prime \prime }(0)}}\right){\displaystyle \frac{1}{\theta ^2}}+{\displaystyle \underset{n=3}{\overset{\mathrm{}}{}}}V(\lambda _n^{(2)})+O(\theta ^3).`$ (30)
In order to estimate the tails of the sums we use the following argument. There is a number $`b(0,s)`$ such that $`V(x)xV^{}(x)`$ for $`0<x<b`$. Taking $`\theta `$ sufficiently large so that $`\lambda _2^{(1)}<b`$ and $`\lambda _3^{(2)}<b`$ we see that

$`{\displaystyle \underset{n=2}{\overset{\mathrm{}}{}}}V(\lambda _n^{(1)})`$ $``$ $`{\displaystyle \underset{n=2}{\overset{\mathrm{}}{}}}\lambda _n^{(1)}V^{}(\lambda _n^{(1)})`$ (31)
$``$ $`\lambda _2^{(1)}{\displaystyle \underset{n=2}{\overset{\mathrm{}}{}}}V^{}(\lambda _n^{(1)})`$ (32)
$``$ $`c\theta ^2{\displaystyle \underset{n=2}{\overset{\mathrm{}}{}}}V^{}(\lambda _n^{(1)}),`$ (33)
where $`c`$ is a constant.

Similarly,

$$\underset{n=3}{\overset{\mathrm{}}{}}V(\lambda _n^{(2)})c\theta ^2\underset{n=3}{\overset{\mathrm{}}{}}V^{}(\lambda _n^{(2)}).$$
(34)
We also have (see )

$$\underset{n=0}{\overset{\mathrm{}}{}}V^{}(\lambda _n^{(1)})=\underset{n=0}{\overset{\mathrm{}}{}}V^{}(\lambda _n^{(2)})=0$$
(35)
so

$$\underset{n=2}{\overset{\mathrm{}}{}}V^{}(\lambda _n^{(1)})=V^{}(\lambda _0^{(1)})V^{}(\lambda _1^{(1)})=O(\theta ^2)$$
(36)
and

$$\underset{n=3}{\overset{\mathrm{}}{}}V^{}(\lambda _n^{(2))}=V^{}(\lambda _0^{(2)})V^{}(\lambda _1^{(2)})V^{}(\lambda _2^{(2)})=O(\theta ^2).$$
(37)
It follows that the tails of the sums in (28) and (30) are $`O(\theta ^4)`$ and Eq. (21) can be written as

$$\left(\frac{1}{V^{\prime \prime }(s)}+\frac{1}{V^{\prime \prime }(0)}\right)\frac{4s^2}{\theta ^2}+O(\theta ^3)=\left(\frac{1}{V^{\prime \prime }(s)}+\frac{1}{V^{\prime \prime }(0)}\right)\frac{8s^2}{\theta ^2}+O(\theta ^3)$$
(38)
which is impossible for large values of $`\theta `$.

By a similar argument one can rule out the existence of a multisoliton solution $`\varphi _{n,z_1,\mathrm{},z_n}`$ at large but finite $`\theta `$ with the properties

$$\varphi _{n,z_1,\mathrm{},z_n}s\underset{i=0}{\overset{n1}{}}|ii|,z_i0$$
(39)
and

$$\varphi _{n,z_1,\mathrm{},z_n}\underset{i=1}{\overset{n}{}}\varphi _{1,z_i},|z_iz_j|1,ij.$$
(40)
The fact that general multisolitons do not exist as solutions to the static equations of motion at finte $`\theta `$ will most likely make it more difficult to establish rigorously the results about moduli space approximation . To actually control such approximations seems to require an analysis of the time dependent equations of motion.

Acknowledgements.  The work of B. D. is supported in part by MaPhySto funded by the Danish National Research Foundation. This research was partly supported by TMR grant no. HPRN-CT-1999-00161. T. J. is indebted to the Niels Bohr Institute for hospitality and would like to thank J. Ambjørn, P. Austing and L. Thorlacius for discussions.


text:
# Hexatic-Herringbone Coupling at the Hexatic Transition in Smectic Liquid Crystals: 4-ϵ Renormalization Group Calculations Revisited

## I Introduction

The nature of phase transitions in two dimensional (2D) systems, has been the subject of numerous investigations over the last three decades. According to the Mermin-Wagner-Hohenberg theorem , the continuous symmetry of the XY and Heisenberg models cannot be spontaneously broken at finite temperature, and there can be no long range magnetic order. However, Kosterlitz and Thouless (KT)  argued that there is a new type of phase transition from a high temperature phase with exponentially decay of the correlations to a low temperature phase with power law decay of the correlations. The idea of KT has been extended by Halperin and Nelson  and Young  to the 2D melting problem. One of the main results of the KTHNY theory is the prediction of an intermediate 2D phase called the hexatic phase for systems that have a six-fold (hexagonal) symmetry in their crystalline ground state. This hexatic phase displays short range positional order, but quasi long range bond-orientational order, which is different from the true long range bond-orientational and quasi long range translational order of a 2D solid phase . The hexatic phase can be characterized by a bond-orientational order parameter defined by $`\mathrm{\Psi }_6=\mathrm{\Psi }_6\mathrm{exp}(i6\psi _6)`$. Assuming that the hexatic state occurs and is not preempted by a direct first order transition from the solid to the isotropic liquid phase, the system should, in the simplest scenario for 2D, display either a KT transition or a first order transition from the hexatic state to the isotropic liquid phase .

It was soon realized after the proposal of the KTHNY theory that the hexatic phases, with short range positional order but true long range bond orientational order, might exist in highly anisotropic three dimensional (3D) systems. Specifically, Birgeneau and Lister  applied the notion of a hexatic state of the 2D melting theory to 3D liquid crystal phases consisting of stacked 2D layers. They proposed that some of the experimentally observed smectic liquid crystal phases could be physical realization of 3D hexatics. Birgeneau and Lister suggested that the (weak) interlayer interaction could promote the quasi long range order of 2D hexatic layers to true long range bond orientational order in 3D.

Stimulated by these theoretical advances, numerous experimental efforts have been undertaken to test theoretical predictions in different liquid crystal materials candidate for displaying hexatic phases . An x-ray study of the liquid crystal compound 65OBC (n-alkyl-4’-m-alkoxybiphenyl-4-carboxylate,n=6,m=5)  provided the first indication of the existence of the 3D analog of the 2D hexatic phase. It was also found that in addition to the hexagonal pattern of diffuse spots of scattered intensity, which is the signature of the hexatic phase, there are some broader peaks corresponding to correlations in the molecular orientations about their long axes . The positions of these peaks show that, locally, the molecules are packed according to a herringbone pattern perpendicular to the smectic layer stacking direction (see Fig. 1). Despite the indication of short range herringbone correlations, this novel phase is simply denoted as the hexatic-B (HexB) phase. Upon increasing temperature, this phase looses its long range bond orientational order and undergoes a transition to the smectic-A (SmA) phase, which essentially consists of a stack of 2D liquid layers. Upon cooling, the HexB phase transforms via a first order phase transition into the crystal-E (CryE) phase, which exhibits both long range translational order and long range herringbone orientational order in the orientations of the molecular axes.

According to the $`U(1)`$ symmetry of the $`\mathrm{\Psi }_6`$ bond orientational order parameter, one would naively expect to find XY-like critical exponents at the SmA$``$HexB transition. However, heat capacity investigations near the SmA$``$HexB transition of 65OBC  and subsequent calorimetric studies on many other components in the nmOBC homologous series  have constantly been reporting continuous (second order) SmA$``$HexB transitions with very large values for the heat capacity critical exponent, $`\alpha 0.6`$. This is drastically different from the 3D XY critical exponent $`\alpha =0.007`$ . As well, thermal conductivity and birefringence experiments have allowed the determination of other static critical exponents, all of which differ systematically from the 3D XY value, while they, together, obey the standard scaling relationships expected for a genuine second order phase transition .

In the light of the existence of the short range herringbone fluctuations, detected in the x-ray diffraction studies , Bruinsma and Aeppli (BA)  formulated a Ginzburg-Landau theory that includes both the hexatic and herringbone order. Because the HexB phase exhibits only short range positional order, BA suggested that the herringbone order can also be represented by an XY order parameter described by $`\mathrm{\Phi }_2=\mathrm{\Phi }_2\mathrm{exp}(i2\varphi _2)`$. At the microscopic level, it is the molecular anisotropy that creates a coupling between the hexatic bond order and the herringbone molecular order . At the phenomenological Ginzburg-Landau level, this coupling is described by a hexatic-herringbone interaction term $`V_{\mathrm{hex}\mathrm{her}}=h\mathrm{Re}(\mathrm{\Psi }_6^{}\mathrm{\Phi }_2^3)`$.

BA constructed an appropriate free energy density based on symmetry considerations and investigated the effects of fluctuation corrections to the mean field behavior for 3D systems. In the mean field approach, their results indicate that the SmA$``$HexB transition should be continuous. However, $`4ϵ`$ renormalization group (RG) calculations, which includes thermal fluctuations and the coupling term $`h\mathrm{Re}(\mathrm{\Psi }_6^{}\mathrm{\Phi }_2^3)`$, show that short-range molecular herringbone correlations coupled to the hexatic ordering drive this transition first order, which becomes second order at a tricritical point.

Interestingly, heat capacity measurement studies of (truly two-dimensional) two-layer free standing films of different nmOBC compounds yield very sharp heat capacity peaks near the SmA$``$HexB transition which can be described by the critical exponent $`\alpha 0.3`$ .This is in sharp contrast with the usual broad and nonsingular specific heat hump associated to the KT transition in the 2D XY model, or yet the first order transition that could occur in a physical system where the vortex core energy is less than some critical value. This $`\alpha 0.3`$ result further suggests that the SmA$``$HexB cannot be described by a simple model with a unique (critical) XY-like order parameter. In this context, there has been some numerical simulations done to obtain more insight into the nature of the SmA$``$HexB transition in 2D systems. The model used in the simulations  consists of a 2D lattice of coupled XY spins based on the BA Hamiltonian. The simulation results suggest the existence of a new type phase transition in which the two different orderings are simultaneously established through a continuous transition. It is interesting to note here that, in seemingly different context, there has also been numerous theoretical and numerical attempts to identify ‘novel chiral’ universality classes for systems such as frustrated XY model and Ising-XY coupled model .

Certainly, for three dimensions, the scenario of a fluctuation-driven first order SmA$``$HexB transition due to hexatic$``$herringbone coupling would appear reasonable for the SmA$``$HexB transition in 65OBC which, upon further cooling, undergoes a HexB$``$CryE transition that establishes long range herringbone and positional order. However, the mixture of 3(10)OBC and PHOAB exhibits a very large temperature range for the HexB phase above the crystallization temperature to the CryE phase. If there were herringbone fluctuations near the SmA$``$HexB transition in that mixture, one could expect them to be quite small because of the large HexB temperature range, and the SmA$``$HexB transition could then possibly be continuous, and to belong to the (then naively expected) XY universality class. However, the fact that the SmA-HexB transition in the 3(10)OBC-PHOAB mixture is found to be first order does not support this simple minded argument . Following the same type of reasoning, recent x-ray diffraction studies on 75OBC  show that the intensity of the herringbone peaks is weaker than those of 65OBC. In principle, if one assumes that 65OBC is near a tricritical point, 75OBC should therefore be further removed from this point (due to the weaker herringbone diffraction peaks, and consequently, weaker $`V_{\mathrm{hex}\mathrm{her}}`$), with again the possibility to recover 3D XY critical behavior. Yet, the same (unconventional) heat-capacity critical exponents are found for these two materials.

The experimental results above could be interpreted as a possible indication of an underlying “novel” (non-XY) stable fixed point that control the SmA$``$HexB transition when herringbone correlations are present. The apparent lack of progress on the theoretical side of the SmA$``$HexB problem has led us to reinvestigate the model of BA and to, specifically, look for a possible calculation error. Firstly, it is important to note that the conclusion of a fluctuation-driven first order within a $`4ϵ`$ calculation is acutely depending on the ‘numerics’ and not constrained by symmetry consideration: a small error (such as a factor 2 instead of 4 here or there) can change the conclusion of a fluctuation-driven first order transition. Secondly, and more specific to the BA problem, we show in the next section, when describing the Ginzburg-Landau free energy density for the SmA$``$HexB transition, that some terms in the RG equations, to first order in $`ϵ`$, were missed in the work of BA. Based on our RG equations, we find two nontrivial strongly-coupled herringbone-hexatic fixed points, apparently missed by those authors. However, those two nontrivial fixed-points are unstable, and we reach the same final conclusion as the one found by BA, namely that of a fluctuation-driven first order transition. We also discuss the possibility of a third and a priori possibly physically pertinent order parameter in the Hamiltonian model of the SmA$``$HexB transition. Because of local distortion of the bond orientational order induced by the anisotropy of the intermolecular potential  and the herringbone correlations, one may generalize the Hamiltonian to the case with three XY-like order parameters, in which two of them are two-fold symmetric, one for the herringbone correlations, $`\mathrm{\Phi }_2`$, and one for the local two-fold distortion, $`\mathrm{\Psi }_2`$, and a third order parameter with six-fold symmetry, $`\mathrm{\Psi }_6`$, related to hexatic ordering. We discuss both mean field and RG calculations for this new model.

The rest of this paper is organized as follows. In Sec. II.A, we reintroduce the BA model and present the result of our RG calculations. In the Sec. II.B, we generalize the Hamiltonian to the case of three order parameter and discuss the mean field theory and RG results. The discussions and conclusions appear in Sec. III.

## II Models and RG Calculations

### A BA Hamiltonian

To formulate the Landau Ginzburg (LG) free energy, which describes both the hexatic and herringbone order, one recalls that the hexatic order is six-fold symmetric, while rotating a herringbone pattern by $`180^o`$ leaves it unchanged. Consequently, the appropriate LG free energy ought to be invariant with respect to the transformation $`\varphi _2(r)\varphi _2(r)+n\pi `$ and $`\psi _6(r)\psi _6(r)+m(2\pi /6)`$ where $`n`$ and $`m`$ are integers. Thus to lowest order in $`\mathrm{\Psi }_6`$ and $`\mathrm{\Phi }_2`$, the BA Hamiltonian  is:

$`\beta F`$ $`=`$ $`{\displaystyle }d^3x[{\displaystyle \frac{r_6}{2}}\mathrm{\Psi }_6^2+{\displaystyle \frac{1}{2}}\mathrm{\Psi }_6^2+{\displaystyle \frac{r_2}{2}}\mathrm{\Phi }_2^2`$ (1)
$`+`$ $`{\displaystyle \frac{1}{2}}\mathrm{\Phi }_2^2+u_6\mathrm{\Psi }_6^4+u_2\mathrm{\Phi }_2^4`$ (2)
$`+`$ $`w\mathrm{\Phi }_2^2\mathrm{\Psi }_6^2+h\mathrm{Re}(\mathrm{\Psi }_6^{}\mathrm{\Phi }_2^3)].`$ (3)
The condition for thermodynamic stability of $`F`$ for $`w=0`$ is $`h^{4/3}<(4^{4/3}/3)u_2u_6^{1/3}`$. This condition can be obtained by minimizing the free energy density on the critical isothermal line $`r_2=r_6=0`$ and requesting that $`\beta F>0`$. As discussed in the Ref. , in the mean field approximation, for $`w=0`$ and $`h=0`$ the phase diagram in the $`r_2r_6`$ plane includes four distinct phases: an isotropic (SmA) phase with $`\mathrm{\Psi }_6=\mathrm{\Phi }_2=0`$, a hexatic (HexB) phase with no herringbone order with $`\mathrm{\Psi }_60,\mathrm{\Phi }_2=0`$, a “putative” ‘herringbone liquid crystal’ phase with $`\mathrm{\Psi }_6=0,\mathrm{\Phi }_20`$, and a fully ordered state with $`\mathrm{\Phi }_20,\mathrm{\Psi }_60`$ akin to a crystalline E phase, with all these phases separated by second order transitions. However, if $`h0`$, the herringbone liquid crystal state with no hexatic order ($`\mathrm{\Phi }_20,\mathrm{\Psi }_6=0`$) is eliminated because $`\mathrm{\Phi }_2`$ acts as a symmetry breaking field on $`\mathrm{\Psi }_6`$. Within mean field theory, the transition lines between the isotropic and ordered phases remain second order for $`h0`$, and terminate together with the first order line at a multicritical point.

We now discuss the RG flow equations and the stability of the fixed points (FPs). Our calculations show that the RG equations to first order in $`ϵ=4d`$ are:

$`{\displaystyle \frac{dr_2}{dl}}=2r_2+{\displaystyle \frac{16K_4u_2}{1+r_2}}+{\displaystyle \frac{4K_4w}{1+r_6}}`$ (4)
$`{\displaystyle \frac{dr_6}{dl}}=2r_6+{\displaystyle \frac{16K_4u_6}{1+r_6}}+{\displaystyle \frac{4K_4w}{1+r_2}}`$ (5)
$`{\displaystyle \frac{du_2}{dl}}=ϵu_240K_4u_2^22K_4w^29K_4h^2`$ (6)
$`{\displaystyle \frac{du_6}{dl}}=ϵu_640K_4u_6^22K_4w^2`$ (7)
$`{\displaystyle \frac{dw}{dl}}=ϵw16K_4wu_216K_4wu_68K_4w^218K_4h^2`$ (8)
$`{\displaystyle \frac{dh}{dl}}=ϵh24K_4hu_212K_4hw,`$ (9)
where $`K_4=1/8\pi ^2`$. The above RG equations differ from those found by BA in Ref.:

* The first set of differences are the $`4K_4w/(1+r_6)`$ and $`4K_4w/(1+r_2)`$ terms in the first and second equations, while BA have $`2K_4w/(1+r_6)`$ and $`2K_4w/(1+r_6)`$. The extra factor $`2`$ comes from the fact that the fields $`\mathrm{\Psi }_6`$ and $`\mathrm{\Phi }_2`$ are complex and the related correlations have two components. For $`h=0`$, Eq. (4) (with the factor $`4K_4w`$) reproduces the RG equations of coupled two components two-vector model as in previous studies. We therefore believe that the above RG equations for $`dr_2/dl`$ and $`dr_6/dl`$ are correct.
* Compared to the BA equations, we also obtain two completely new terms, $`18K_4h^2`$, in the fifth equation and $`12K_4hw`$ in the sixth equation, which can be easily checked using Feynmann diagram technique. Indeed, these two new terms come from the connected diagrams in the second order perturbative RG obtained by multiplication of the relevant diagrams of $`h\mathrm{\Psi }_6^{}\mathrm{\Phi }_2^3`$ with $`h\mathrm{\Psi }_6\mathrm{\Phi }_2^3`$, and of $`w\mathrm{\Phi }_2^2\mathrm{\Psi }_6^2`$ with $`h\mathrm{Re}(\mathrm{\Psi }_6^{}\mathrm{\Phi }_2^3)`$, respectively.

Because of the two extra terms in the RG equations for $`dw/dl`$ and $`dh/dl`$, we obtain, in addition to the simple decoupled FP ($`r_6^{}=r_2^{}=ϵ/5,u_6^{}=u_2^{}=ϵ/(40K_4)`$, $`w^{}=h^{}=0`$), two new fixed points such that $`(w^{}0,h^{}=0)`$ and ($`w^{}0,h^{}0)`$. The first nontrivial FP is given by $`h^{}=0,r_6^{}=r_2^{}=ϵ/4,u_6^{}=u_2^{}=ϵ/(48K_4)`$, and w$`{}_{}{}^{}=ϵ/(24K_4)`$. This FP, akin to the one found in minimally coupled two component two vector model , was not discussed by BA.

However, and most interestingly, we find another “new” nontrivial mixed herringbone-hexatic FP with all the couplings being non-zero:

$`r_6^{}`$ $`=`$ $`0.24845566ϵ,`$ (10)
$`r_2^{}`$ $`=`$ $`0.24018995ϵ,`$ (11)
$`u_6^{}`$ $`=`$ $`0.01941403ϵ/K_4,`$ (12)
$`u_2^{}`$ $`=`$ $`0.01838082ϵ/K_4,`$ (13)
$`w^{}`$ $`=`$ $`0.04657169ϵ/K_4,`$ (14)
$`h^{}`$ $`=`$ $`\pm 0.00766519ϵ/K_4.`$ (15)
Therefore, based on our RG calculations, there is a FP with $`h^{}0`$, which was not found in the previous work of BA. Linearizing the recursion relations in the vicinity of the FPs yields for the FP with ($`w^{}0`$, $`h^{}=0`$): $`y_1=2ϵ/2,y_2=2ϵ/6,y_3=ϵ,y_4=2ϵ/3,y_5=y_6=0`$. Thus, two eigenvalues are marginal, compatible with what has been found in similar minimally coupled two component two vector model . The eigenvalues for the FP with ($`h^{}0,w^{}0)`$ above are

$`y_1`$ $`=`$ $`20.488829ϵ,`$ (17)
$`y_2`$ $`=`$ $`20.115889ϵ,`$ (18)
$`y_3`$ $`=`$ $`0.997894ϵ,`$ (19)
$`y_4`$ $`=`$ $`0.537266ϵ,`$ (20)
$`y_5`$ $`=`$ $`+0.121467ϵ,`$ (21)
$`y_6`$ $`=`$ $`+0.0402392ϵ.`$ (22)
These results show that there are four positive eigenvalues, and the above “new” nontrivial FP is therefore unstable. The two largest (most positive) eigenvalues, $`y_1`$ and $`y_2`$, correspond, respectively, to the thermal eigenvalue and the ‘relative’ coupling strength that positions the system in coupling parameter space and determines what sequence of phase transition occurs. Namely, isotropic$``$(hexatic+herringbone) via a unique phase transtion or isotropic $``$ hexatic $``$(hexatic+herringbone) via two distinct phase transitions. The four eigenvalues $`y_3`$$``$$`y_6`$ essentially control the flow in the $`wh`$ plane. $`y_5`$ and $`y_6`$ are positive, and we interpret the result as an indication that the above “new” nontrivial mixed herringbone-hexatic FP is unstable. Indeed, we have confirmed this by explicit numerical integration of the RG equations, and found that the RG flow goes to the unstable region identified above, which we interpret as the transition being driven first order by fluctuations. Therefore, while we have indeed found some discrepancies between our RG equations and those of BA, and recovered two extra coupled fixed points, we at the end still reach the same physical conclusion of BA, namely that of a fluctuation-driven first order smA$``$HexB transition.

### B Generalized Hamiltonian

We expect physically that the local molecular anisotropy (e.g. from the anisotropic nature of benzene rings found in most thermotropic liquid crystal materials) to couple to the local bond direction, and to create a local two-fold distortion of the otherwise perfect local six-fold symmetric nearest-neighbor bond order. Consequently, we now discuss the LG free energy, which describes both the hexatic and herringbone order, as well as the local two-fold distortion of the bond order. If we assume that the distortion of lattice has two-fold symmetry, with the order parameter $`\mathrm{\Psi }_2=\mathrm{\Psi }_2\mathrm{exp}(i2\psi _2)`$, then the resulting free energy is invariant under the transformation $`\varphi _2(r)\varphi _2(r)+n\pi `$, $`\psi _6(r)\psi _6(r)+m(2\pi /6)`$, and $`\psi _2(r)\psi _2(r)+p\pi `$ where $`n`$, $`m`$, and $`p`$ are integers.

From an RG point of view, our motivation to expand the symmetry of our Hamiltonian stems from the observation that in $`N`$ coupled 2-vector models a new stable fixed points (called mixed fixed point ) appears in the coupling parameter space (when $`N>2`$). Thus, to lowest order in $`\mathrm{\Psi }_6`$, $`\mathrm{\Phi }_2`$, and $`\mathrm{\Psi }_2`$, we have $`\beta F=\beta F_0+U`$, where the Gaussian part is given by

$`\beta F_0`$ $`=`$ $`{\displaystyle \frac{1}{2}}{\displaystyle }d^3x[r_6\mathrm{\Psi }_6^2+r_2\mathrm{\Phi }_2^2+\stackrel{~}{r}_2\mathrm{\Psi }_2^2`$ (23)
$`+`$ $`2r\mathrm{Re}(\mathrm{\Phi }_2\mathrm{\Psi }_2^{})+\mathrm{\Psi }_6^2+\mathrm{\Phi }_2^2`$ (24)
$`+`$ $`\mathrm{\Psi }_2^2+2g\mathrm{Re}(\mathrm{\Phi }_2\mathrm{\Psi }_2^{})],`$ (25)
and the perturbative Hamiltonian has the following form:

$`U`$ $`=`$ $`{\displaystyle }d^3x[u_6\mathrm{\Psi }_6^4+u_2\mathrm{\Phi }_2^4+\stackrel{~}{u}_2\mathrm{\Psi }_2^4`$ (26)
$`+`$ $`w_1\mathrm{\Psi }_6^2\mathrm{\Phi }_2^2+w_2\mathrm{\Psi }_6^2\mathrm{\Psi }_2^2+w_3\mathrm{\Phi }_2^2\mathrm{\Psi }_2^2`$ (27)
$`+`$ $`h_1\mathrm{Re}(\mathrm{\Psi }_6^{}\mathrm{\Phi }_2^3)+h_2\mathrm{Re}(\mathrm{\Psi }_6^{}\mathrm{\Psi }_2^3)+h_3\mathrm{Re}(\mathrm{\Phi }_2^2\mathrm{\Psi }_2^2)`$ (28)
$`+`$ $`v_1\mathrm{Re}(\mathrm{\Psi }_6^{}\mathrm{\Phi }_2\mathrm{\Psi }_2^2)+v_2\mathrm{Re}(\mathrm{\Psi }_6^{}\mathrm{\Phi }_2^2\mathrm{\Psi }_2)+v_3\mathrm{\Psi }_6^2\mathrm{Re}(\mathrm{\Phi }_2\mathrm{\Psi }_2^{})`$ (29)
$`+`$ $`v_4\mathrm{\Phi }_2^2\mathrm{Re}(\mathrm{\Phi }_2\mathrm{\Psi }_2^{})+v_5\mathrm{\Psi }_2^2\mathrm{Re}(\mathrm{\Phi }_2\mathrm{\Psi }_2^{})]`$ (30)
For the case that $`r_2=\stackrel{~}{r}_2`$, one can simply diagonalize the Gaussian part of the Hamiltonian using the transformation $`\mathrm{\Phi }_2=(\stackrel{~}{\mathrm{\Phi }}_2+\stackrel{~}{\mathrm{\Psi }}_2)/\sqrt{2}`$ and $`\mathrm{\Psi }_2=(\stackrel{~}{\mathrm{\Phi }}_2\stackrel{~}{\mathrm{\Psi }}_2)/\sqrt{2}`$, and then do the RG calculations. The RG calculations for the case that $`\mathrm{\Psi }_6=0`$ were done by Yosefin and Domany  in the study of the phase transitions in fully frustrated XY models

In the mean field approach, there are now four distinct phases for $`r0`$ ($`\mathrm{\Psi }_6=\mathrm{\Phi }_2=\mathrm{\Psi }_2=0;\mathrm{\Psi }_6=0,\mathrm{\Phi }_20,\mathrm{\Psi }_20;\mathrm{\Psi }_60,\mathrm{\Psi }_2=\mathrm{\Phi }_2=0;\mathrm{\Psi }_60,\mathrm{\Psi }_20,\mathrm{\Phi }_20`$). One should note that for the phases where both $`\mathrm{\Phi }_2`$ and $`\mathrm{\Psi }_2`$ are non-zero, we have the condition of local stability of the free energy, $`r_2\stackrel{~}{r}_2<r^2`$. In addition, the singularity of the propagators at zero wavevector ($`𝐪)=0`$) is for $`r_2\stackrel{~}{r}_2=r^2`$, which is the critical point of the system.

To obtain further insight into the specific situation where both $`\mathrm{\Phi }_2`$ and $`\mathrm{\Psi }_2`$ go simultaneously critical (soft), we perform an RG calculation. To simplify the calculations, we consider the case that $`r_2=\stackrel{~}{r}_2`$, so that the fields $`\mathrm{\Phi }_2`$ and $`\mathrm{\Psi }_2`$ are simultaneously critical (soft), and that they are equally coupled to the $`\mathrm{\Psi }_6`$ hexatic field ($`u_2=\stackrel{~}{u}_2`$, $`w_1=w_2`$). We further need to require that the full theory, after diagonalization of the gaussian part, is self-consistent with no new RG-generated terms. This imposes that $`h_1=h_2`$ and $`v_1=v_2`$. Using the above mentioned transformation for $`\mathrm{\Phi }_2`$ and $`\mathrm{\Psi }_2`$ and rescaling the fields, one can rewrite the Hamiltonian as,

$`\beta F={\displaystyle }d^3x[{\displaystyle \frac{r_6}{2}}\mathrm{\Psi }_6^2+{\displaystyle \frac{r_2^{}}{2}}\stackrel{~}{\mathrm{\Phi }}_2^2+{\displaystyle \frac{\stackrel{~}{r}_2^{}}{2}}\stackrel{~}{\mathrm{\Psi }}_2^2`$ (31)
$`+{\displaystyle \frac{1}{2}}\mathrm{\Phi }_6^2+{\displaystyle \frac{1}{2}}\stackrel{~}{\mathrm{\Phi }}_2^2+{\displaystyle \frac{1}{2}}\stackrel{~}{\mathrm{\Psi }}_2^2+u_6\mathrm{\Psi }_6^4`$ (32)
$`+u_2^{}\stackrel{~}{\mathrm{\Phi }}_2^4+\stackrel{~}{u}_2^{}\stackrel{~}{\mathrm{\Psi }}_2^4+w_1^{}\mathrm{\Psi }_6^2\stackrel{~}{\mathrm{\Phi }}_2^2`$ (33)
$`+w_2^{}\mathrm{\Psi }_6^2\stackrel{~}{\mathrm{\Psi }}_2^2+w_3^{}\stackrel{~}{\mathrm{\Phi }}_2^2\stackrel{~}{\mathrm{\Psi }}_2^2`$ (34)
$`+h_1^{}\mathrm{Re}(\mathrm{\Psi }_6^{}\stackrel{~}{\mathrm{\Phi }}_2^3)+h_3^{}\mathrm{Re}(\stackrel{~}{\mathrm{\Phi }}_2^2\stackrel{~}{\mathrm{\Psi }}_2^2)+v_1^{}\mathrm{Re}(\mathrm{\Psi }_6^{}\stackrel{~}{\mathrm{\Phi }}_2\stackrel{~}{\mathrm{\Psi }}_2^2)],`$ (35)
where the new (primed) coefficients can be written in terms of old coefficients. To first order of $`ϵ`$, the RG equations are given by,

$`{\displaystyle \frac{dr_6}{dl}}`$ $`=`$ $`2r_6+{\displaystyle \frac{16K_4u_6}{1+r_6}}+{\displaystyle \frac{4K_4w_1^{}}{1+r_2^{}}}+{\displaystyle \frac{4K_4w_2^{}}{1+\stackrel{~}{r}_2^{}}}`$ (36)
$`{\displaystyle \frac{dr_2^{}}{dl}}`$ $`=`$ $`2r_2^{}+{\displaystyle \frac{16K_4u_2^{}}{1+r_2^{}}}+{\displaystyle \frac{4K_4w_1^{}}{1+r_6}}+{\displaystyle \frac{4K_4w_3^{}}{1+\stackrel{~}{r}_2^{}}}`$ (37)
$`{\displaystyle \frac{d\stackrel{~}{r}_2^{}}{dl}}`$ $`=`$ $`2\stackrel{~}{r}_2^{}+{\displaystyle \frac{16K_4\stackrel{~}{u}_2^{}}{1+\stackrel{~}{r}_2^{}}}+{\displaystyle \frac{4K_4w_2^{}}{1+r_6}}+{\displaystyle \frac{4K_4w_3^{}}{1+r_2^{}}}`$ (38)
$`{\displaystyle \frac{du_6}{dl}}`$ $`=`$ $`ϵu_640K_4u_6^22K_4w_1^22K_4w_2^2`$ (39)
$`{\displaystyle \frac{du_2^{}}{dl}}`$ $`=`$ $`ϵu_2^{}40K_4u_2^22K_4w_1^22K_4w_3^29K_4h_1^2`$ (40)
$``$ $`2K_4h_3^2`$ (41)
$`{\displaystyle \frac{d\stackrel{~}{u}_2^{}}{dl}}`$ $`=`$ $`ϵ\stackrel{~}{u}_2^{}40K_4\stackrel{~}{u}_2^22K_4w_2^22K_4w_3^22K_4h_3^2`$ (42)
$``$ $`K_4v_1^2`$ (43)
$`{\displaystyle \frac{dw_1^{}}{dl}}`$ $`=`$ $`ϵw_1^{}16K_4u_6w_1^{}16K_4u_2^{}w_1^{}8K_4w_1^2`$ (44)
$``$ $`4K_4w_2^{}w_3^{}18K_4h_1^22K_4v_1^2`$ (45)
$`{\displaystyle \frac{dw_2^{}}{dl}}`$ $`=`$ $`ϵw_2^{}16K_4u_6w_2^{}16K_4\stackrel{~}{u}_2^{}w_2^{}8K_4w_2^2`$ (46)
$``$ $`4K_4w_1^{}w_3^{}4K_4v_1^2`$ (47)
$`{\displaystyle \frac{dw_3^{}}{dl}}`$ $`=`$ $`ϵw_3^{}16K_4u_2^{}w_3^{}16K_4\stackrel{~}{u}_2^{}w_3^{}8K_4w_3^2`$ (48)
$``$ $`4K_4w_1^{}w_2^{}16K_4(h_3^{})^24K_4v_1^2`$ (49)
$`{\displaystyle \frac{dh_1^{}}{dl}}`$ $`=`$ $`ϵh_1^{}24K_4h_1^{}u_2^{}12K_4w_1^{}h_1^{}4K_4v_1^{}h_3^{}`$ (50)
$`{\displaystyle \frac{dh_3^{}}{dl}}`$ $`=`$ $`ϵh_3^{}8K_4u_2^{}h_3^{}8K_4\stackrel{~}{u}_2^{}h_3^{}16K_4w_3^{}h_3^{}`$ (51)
$``$ $`6K_4h_1^{}v_1^{}`$ (52)
$`{\displaystyle \frac{dv_1^{}}{dl}}`$ $`=`$ $`ϵv_1^{}8K_4v_1^{}\stackrel{~}{u}_2^{}4K_4w_1^{}v_1^{}8K_4w_2^{}v_1^{}`$ (53)
$``$ $`8K_4w_3^{}v_1^{}12K_4h_1^{}h_3^{}.`$ (54)
We have found a number of FPs for the above RG equations that fulfill the condition ($`r_2^{}=\stackrel{~}{r}_2^{}`$, and for which $`u_2^{}=\stackrel{~}{u}_2^{}`$, $`w_1^{}=w_2^{}`$, $`h_1^{}=h_2^{}`$, $`v_1^{}=v_2^{}`$ and $`v_4^{}=v_5^{}`$). Some of these are given in the Table I, where each column correspond to a different FP (the nontrivial mixed herringbone-hexatic FP found in Section IIa occurs here as well, but is again unstable). As in the simplest case of Hamiltonian (2.1), none of the FP correspond to a nontrivial stable fixed point: each fixed point (column) shows more than two positive eigenvalues. Given the complexity of those nonlinear equations we cannot be 100% sure that we have found all the (unstable) FP in the $`(u_6,u_6^{},\stackrel{~}{u}_2^{},w_1^{},w_3^{},h_1^{},h_3^{},v_1^{})`$ plane. However, a numerical investigation of the RG flow in that plane starting from a large number of initial values for the coupling parameters always failed to converge towards a stable (attractive) fixed point. Therefore, we reach the same conclusion obtained in the previous section, namely that of a fluctuation-driven first order transition for this expanded symmetry Hamiltonian, and the stabilization of novel fixed point in $`N`$ coupled two vector models for $`N>2`$ does not appear to occur in this generalized Hamiltonian due to the extra coupling parameters.

Table I: The FPs of RG equations for generalized Hamiltonian.

| $`r_6^{}/ϵ`$ | -1/5 | -3/11 | -2/7 | 0 | -0.2761944898 | -0.2583527792 | -0.2727272730 |
| --- | --- | --- | --- | --- | --- | --- | --- |
| $`r_2^{}/ϵ`$ | -1/5 | -3/11 | -2/7 | -1/5 | -0.2761944898 | -0.2756996328 | -0.2727272730 |
| $`\stackrel{~}{r}_2^{}/ϵ`$ | -1/5 | -3/11 | -2/7 | -1/5 | -0.2566490108 | -0.2756996328 | -0.2727272724 |
| $`u_6^{}(K_4/ϵ)`$ | 1/40 | 1/44 | 1/56 | 0 | 0.0202792491 | 0.0237809941 | 0.0227272727 |
| $`u_2^{}(K_4/ϵ)`$ | 1/40 | 1/44 | 1/56 | 0.0125 | 0.0202792491 | 0.0194412033 | 0.0170454546 |
| $`\stackrel{~}{u}_2^{}(K_4/ϵ)`$ | 1/40 | 1/44 | 1/56 | 0.0125 | 0.0238702512 | 0.0194412033 | 0.0170454546 |
| $`w_1^{}(K_4/ϵ)`$ | 0 | 1/44 | 1/28 | 0 | 0.0405584982 | 0.0170262066 | 0.0227272727 |
| $`w_2^{}(K_4/ϵ)`$ | 0 | 1/44 | 1/28 | 0 | 0.0164217503 | 0.0170262066 | 0.0227272727 |
| $`w_3^{}(K_4/ϵ)`$ | 0 | 1/44 | 1/28 | 0.05 | 0.0164217503 | 0.0430587966 | 0.0454545454 |
| $`h_1^{}(K_4/ϵ)`$ | 0 | 1/44 | 0 | 0 | 0 | 0 | 0 |
| $`h_3^{}(K_4/ϵ)`$ | 0 | 1/44 | 0 | $`\pm `$ 0.025 | 0 | $`\pm `$ 0.0041763898 | $`\pm `$ 0.0113636364 |
| $`v_1^{}(K_4/ϵ)`$ | 0 | 1/44 | 0 | 0 | 0 | 0 | 0 |

## III Discussion

In Section II. A, we found that there are two “new” FPs missed in Ref. . However, those FPs are unstable, and we reach the same final conclusion as Bruinsma and Aeppli, namely that the SmA$``$HexB transition is driven first order by fluctuations in the BA Hamiltonian. We also discussed a slightly modified simple model that considers not only the hexatic and herringbone order, but also one that involves the local two-fold deformation of the bond correlations induced by the herringbone correlations. We assumed that this deformation can be represented by a two-fold symmetry order parameter, as in the case of herringbone order, and wrote the free energy based on symmetry arguments. We were not able to find a new stable FP which could possibly result in unconventional (“new”) second order critical exponents.

It is not clear to us in what directions to pursue the paradoxical puzzle of “new universality” (non $`3D`$ XY transition) in hexatic liquid crystal materials. In the work presented here, we have found a novel “mixed” hexatic-herringbone fixed point in the theory, but which is unstable to first order in $`ϵ`$. One note that the two positives eigenvalues $`y_5`$ and $`y_6`$ in Eq. 2.4 are only very slightly positive for $`ϵ=1`$. This observation may open the possibility that in a calculation that includes hexatic-herringbone coupling, there is no stable fixed point to lowest order in $`ϵ`$, but that the fixed point may be stabilized in a theory that goes beyond an $`O(ϵ^1)`$ calculation. This is what happens, for example, in the normal to superconducting phase transition where the $`ϵ`$ expansion to order $`ϵ^1`$ predicts a fluctuation-driven first order phase transition while theoretical arguments and Monte Carlo simulations strongly argue for a second order (inverted 3D XY) phase transition.

.

Having said that, one should note that there are experiments on liquid crystal materials that do not display any herringbone correlations  but still show a SmA$``$HexB transition with exponents that differ from the critical behavior expected to a $`2D`$  or $`3D`$ XY critical behavior. That may suggest that the whole idea of hexatic$``$herringbone interactions is a red herring. Another possibility is that of (more subtle) ‘hidden’ order parameter(s) distinct from the herringbone order characterizes the SmA$``$HexB transition in real materials, and that the coupling between this hidden order parameter and the hexatic order parameter $`\mathrm{\Psi }_6`$ produces a novel stable fixed point. Clearly, more experimental studies are needed to shed light on this problem. In particular, high resolution scattering experiments would seem necessary to search for extended short-range correlation in either molecular correlations or distortion of hexagonal co-ordination to shed some light on what such ‘hidden’ order parameter(s) may be.

We finally note that a related (unconventional critical behavior) situation arises in the context of layered systems of smectic liquid crystals studied by Defontaines and Prost . These authors have argued that critical points that do not involve any symmetry change can define a set of new universality classes in layered systems. Possibly, considerations of some features of the Defontaines and Prost theory may be useful in further investigations of the SmA$``$HexB problem.

We hope that our work and reinvestigation of the long standing SmA$``$Hex-B transition in smectic liquid crystal materials will motivate further theoretical, numerical and experimental investigations of this very interesting problem.

## IV Acknowledgements

We thank Ian Affleck, Igor Herbut and and Mehran Kardar for useful discussions. M. G. thanks C. C. Huang for stimulating his interest in the problem of the critical behavior of the SmA$``$HexB transition. This research has been funded by NSERC of Canada. M. G. acknowledges the Research Corporation and the Province of Ontario for generous financial support.


text:
# Dynamical Screening and Superconducting State in Intercalated Layered Metallochloronitrides

## Abstract

An essential property of layered systems is the dynamical nature of the screened Coulomb interaction. Low energy collective modes appear as a consequence of the layering and provide for a superconducting-pairing channel in addition to the electron-phonon induced attractive interaction. We show that taking into account this feature allows to explain the high critical temperatures ($`T_c26`$K) observed in recently discovered intercalated metallochloronitrides. The exchange of acoustic plasmons between carriers leads to a significant enhancement of the superconducting critical temperature that is in agreement with the experimental observations.

preprint: Phys. Rev. B 66, 100501(R) (2002)

Screening of the Coulomb interaction takes very different forms in layered conductors and three dimensional (3D) isotropic metals. We show that the dynamic screening in layered systems can lead to a Coulomb induced enhancement of the superconducting pairing and might be an essential addition to the usual electron-phonon contribution. This important feature results from the existence of low energy electronic collective modes characteristic for layered materials.

The aim of the present paper is to explain the nature of the superconducting state in layered intercalated metallochloronitrides yamanaka . It has been shown that intercalation of metallic ions and organic molecules into the parent compound (Zr,Hf)NCl leads to a superconductor with rather high critical temperature ($`T_c26`$K) yamanaka . Based on experimental studies yamanaka ; kawaji ; shamoto1 ; shamoto2 ; adelmann ; tou1 ; yokoya  and band structure calculations weht ; hase  it was concluded that 1) electron-phonon mediated pairing is insufficient to explain the high $`T_c`$s observed and that 2) there is no evidence for the presence of strong correlations; the system can be described within Fermi liquid theory. In addition, these compounds do not have magnetic ions which excludes a magnetic mechanism as well. No explanation has been suggested sofar as to what pairing mechanism can allow to reach the observed critical temperatures. The theory proposed below shows that such high $`T_c`$s can be obtained by including the additional pairing contribution arising from the interaction of carriers with acoustic plasmons; this is the manifestation of the dynamic screening effect of the Coulomb interaction.

The description of layered conductors can be made by neglecting the small interlayer hopping in first approximation. On the other hand, it is essential to take into account the screened interlayer Coulomb interaction which has an important dynamic part. Indeed, it is known that for usual 3D materials this interaction can be considered in the static limit since electronic collective modes are very high in energy (the optical plasmon energies are of the order 5-30eV in metals; see,e.g., nozieres .) Therefore, the Coulomb repulsion enters the theory of superconductivity as a single constant pseudopotential $`\mu ^{}`$. The situation is very different in layered conductors: incomplete screening of the Coulomb interaction results from the layering fetter . The response to a charge fluctuation is time-dependent and the frequency dependence of the screened Coulomb interaction becomes important. This leads to the presence of low energy electronic collective modes: the acoustic plasmons. It is this particular feature of layered materials that brings about an additional contribution to the pairing interaction between electrons.

The order parameter $`\mathrm{\Delta }(𝐤,\omega _n)`$ of the superconducting state is described by

$`\mathrm{\Delta }(𝐤,\omega _n)Z(𝐤,\omega _n)`$
$`=T{\displaystyle \underset{m=\mathrm{}}{\overset{\mathrm{}}{}}}{\displaystyle \frac{\mathrm{d}^3𝐤^{}}{(2\pi )^3}\mathrm{\Gamma }(𝐤,𝐤^{};\omega _n\omega _m)F^{}(𝐤,\omega _m)},`$
where $`F^{}=<c_{k,}^{}c_{k,}^{}>`$ is the Gor’kov pairing function) and $`Z(k,\omega _n)`$ is the renormalization function , defined by

$`Z(𝐤,\omega _n)1`$
$`={\displaystyle \frac{T}{\omega _n}}{\displaystyle \underset{m=\mathrm{}}{\overset{\mathrm{}}{}}}{\displaystyle \frac{\mathrm{d}^3𝐤^{}}{(2\pi )^3}\mathrm{\Gamma }(𝐤,𝐤^{};\omega _n\omega _m)G(𝐤,\omega _m)}.`$
$`G=<c_{k,\sigma }^{}c_{k,\sigma }>`$ is the usual Green function, and $`\mathrm{\Gamma }`$ the total interaction kernel; $`\omega _n=(2n+1)\pi T`$. We use the thermodynamic Green’s functions formalism (see e.g. Ref. AGD ). The $`T_c`$ for layered superconductors is obtained by solving the set of equations (Dynamical Screening and Superconducting State in Intercalated Layered Metallochloronitrides,Dynamical Screening and Superconducting State in Intercalated Layered Metallochloronitrides) self-consistently.

The interaction kernel is composed of two parts, $`\mathrm{\Gamma }=\mathrm{\Gamma }_{ph}+\mathrm{\Gamma }_c`$, where

$`\mathrm{\Gamma }_{ph}(𝐪;|nm|)`$ $`=`$ $`|g_\nu (𝐪)|^2D(𝐪,|nm|)`$ (3)
$`\mathrm{\Gamma }_c(𝐪;|nm|)`$ $`=`$ $`{\displaystyle \frac{V_c(𝐪)}{ϵ(𝐪,|nm|)}}.`$ (4)
The first term, $`\mathrm{\Gamma }_{ph}`$, is the usual pairing contribution resulting from the electron-phonon interaction. $`D(q,nm)=\mathrm{\Omega }_\nu ^2(q)[\omega _n\omega _m)^2+\mathrm{\Omega }_\nu ^2(q)]^1`$ is the phonon temperature Green’s function, $`\mathrm{\Omega }_\nu (q)`$ the phonon frequency; summation over phonon branches $`\nu `$ is assumed. The second contribution to the interaction kernel, $`\mathrm{\Gamma }_c`$, is the Coulomb part written in its most general form as the ratio of the bare Coulomb interaction $`V_c(𝐪)`$ and the dielectric function $`ϵ(𝐪,\omega _n\omega _m)`$. Both functions have to be calculated for a layered structure.

The Coulomb interaction for conducting layers separated by spacers of dielectric constant $`ϵ_M`$ can be written in the form fetter ; kresin :

$`V_c(𝐪)={\displaystyle \frac{2\pi e^2}{ϵ_Mq_{}}}R(q_{},q_z)={\displaystyle \frac{\lambda _c}{N(E_F)}}\stackrel{~}{V}_c(q_{};q_z),`$ (5)
where $`q_{}(q_z)`$ is the in-plane (out-of-plane) component of the wave-vector. In the last expression, $`N(E_F)`$ is the 2D electronic density of states (DoS) at the Fermi energy $`E_F`$, and $`\stackrel{~}{V}_c(𝐪_{};q_z)=R(q_{},q_z)/(2k_Fq_{})`$ with

$`R(𝐪_{},q_z)={\displaystyle \frac{\mathrm{sinh}(q_{}L)}{\mathrm{cosh}(q_{}L)\mathrm{cos}(q_zL)}}.`$ (6)
$`L`$ is the interlayer spacing and $`\lambda _c=(e^2/\mathrm{}v_F)/2ϵ_M`$ is the dimensionless Coulomb interaction constant ($`v_F`$ is the Fermi velocity). The dielectric function $`ϵ(𝐪,\omega _n\omega _m)`$ has been calculated for a layered system kresin ; morawitz ; bill99  in RPA. It has been shown there that the plasmon spectrum contains anisotropic bands $`\omega _{pl}=\omega _{pl}(𝐪_{},q_z)`$ that can be labeled by $`q_z`$ and which are the low frequency acoustic modes.

Eqs. (Dynamical Screening and Superconducting State in Intercalated Layered Metallochloronitrides,Dynamical Screening and Superconducting State in Intercalated Layered Metallochloronitrides) can be cast into the following matrix form near $`T_c`$ (see our paper, Ref. bill99 )

$`{\displaystyle \underset{m}{}}{\displaystyle \underset{k_z^{}}{}}K_{n,m}(|q_z|)\mathrm{\Phi }_m(k_z^{})=\eta \mathrm{\Phi }_n(k_z),`$ (7)
where $`q_zk_zk_z^{}`$ are the wave-vector components normal to the conducting layers and $`\mathrm{\Phi }_m(k_z^{})=\mathrm{\Delta }_m(k_z^{})/\sqrt{2m+1}`$ is the reduced order parameter. In the case of a layered superconductor, the matrix K takes the form

$`K_{n,m}(|q_z|)`$ $`=`$ $`{\displaystyle \frac{1}{N_z}}{\displaystyle \frac{1}{\sqrt{2n+1}\sqrt{2m+1}}}\{`$
$`\lambda \text{}\left[D(nm)+D(n+m+1)\right]`$
$`+\lambda _c\left[\mathrm{\Gamma }_c^I(|nm|;|q_z|)+\mathrm{\Gamma }_c^I(|n+m+1|;|q_z|)\right]`$
$`\mu ^{}\theta (\mathrm{\Omega }_c|\omega _m|)`$
$`\delta _{n,m}{\displaystyle \underset{p=0}{\overset{2n}{}}}[\lambda D(np)+\lambda _c\mathrm{\Gamma }_c^I(|np|;|q_z|)]\}.`$
$`nm`$ is short hand for the difference of Matsubara frequencies $`\omega _n\omega _m=2\pi \stackrel{~}{T}(nm)`$ \[with $`\stackrel{~}{T}=k_BT/\mathrm{\Omega }`$; we consider an Einstein phonon $`\mathrm{\Omega }_\nu (𝐪)\mathrm{\Omega }`$\]. $`\mathrm{\Omega }_c`$ is the cutoff used to define the pseudopotential $`\mu ^{}`$, and $`N_z`$ is the number of $`q_z`$-points considered in the Brillouin zone. All but the static Coulomb repulsion $`\mu ^{}`$ are temperature dependent quantities. The critical temperature of the superconducting phase-transition $`T_c`$ is reached when the highest eigenvalue is $`\eta =1`$.

$`\mathrm{\Gamma }_c^I(n,|q_z|)`$ is the frequency-dependent contribution of the screened Coulomb interaction arising from acoustic plasmons. It has been shown by two of the authors, I. Bozovic, G. Rietveld, and D. van der Marel morawitz  that the density of states (DoS) of the low-energy collective modes is peaked at $`q_z=\pi `$ and $`q_z=0`$. Furthermore, it was demonstrated that the $`q_z=0`$ term is repulsive and can therefore be included into the pseudopotential $`\mu ^{}`$ bill99 . The main plasmon contribution to the pairing is thus obtained for $`q_z=\pi /L`$ and has the form

$`\mathrm{\Gamma }_c^I(nm)`$ $`=`$ $`{\displaystyle \frac{\lambda _c}{2\pi }}{\displaystyle _0^{\stackrel{~}{q}_c}}{\displaystyle \frac{\mathrm{d}\stackrel{~}{\mathrm{q}}}{\sqrt{1\stackrel{~}{\mathrm{q}}^2}}}{\displaystyle \frac{\stackrel{~}{V}_c(\stackrel{~}{q})}{ϵ(\stackrel{~}{q},nm)}},`$ (9)
where $`\stackrel{~}{q}q_{}/2k_F`$, and $`\stackrel{~}{q}_c=\mathrm{min}\{1,|\omega _n\omega _m|/4E_F\}`$ divides $`(\omega ,q)`$-space into the regions $`\omega >qv_F`$ and $`\omega <qv_F`$. The first region corresponds to the dynamic response and contains plasmon excitations, including the acoustic plasmon branches. In the second region the response can be treated in the static approximation and represents the usual repulsive part of the screened Coulomb interaction. We calculate the value of the critical temperature from Eqs. (7-9).

In order to demonstrate the importance of dynamic screening for superconductivity we calculate $`T_c`$ for the following set of realistic parameters: $`L=15`$Å, $`\lambda =0.5`$, $`\mathrm{\Omega }=70`$meV, $`ϵ_M=3`$, $`v_F=5\times 10^7`$, $`\mu ^{}=0.1`$, $`m^{}=m_e`$. As will be seen below, these values are close to those found in metallochloronitrides. With these parameters, the Coulomb interaction constant defined earlier is $`\lambda _c=0.6`$. One can, therefore, use RPA in a first approximation and neglect vertex corrections.

With use of aforementioned values for the three quantities $`\lambda `$, $`\mathrm{\Omega }`$ and $`\mu ^{}`$ one can, in a first step, calculate the value $`T_{c,ph}`$ which is the critical temperature in the absence of dynamic screening ($`\mathrm{\Gamma }_c^I=0`$). One obtains $`T_{c;ph}=12`$K. If we now take into account the effect of dynamic screening and calculate $`T_c`$ using all parameters given above, we obtain $`T_c=25`$K. This demonstrates that the value of $`T_c`$ in layered superconductors can be drastically affected (enhanced) by the dynamic part of the screened Coulomb interaction.

We now apply our approach to a specific case among intercalated metallochloronitrides. Namely the compound Li0.48(THF)yHfNCl which has a $`T_c=25.5`$K yamanaka . We selected this compound as a study case because there has been relatively detailed experimental and theoretical work done on this layered material. From Ref. adelmann ; tou1  the interlayer distance $`L`$ and characteristic phonon frequency $`\mathrm{\Omega }`$ are equal to $`L=18.7`$Å and $`\mathrm{\Omega }=60`$meV, respectively. The effective mass and Fermi energy have been estimated from band structure calculations weht . Accordingly, $`m^{}/m_e0.6`$ where $`m_e`$ is the free electron mass and $`E_F1`$eV. Finally, according to Ref. tou1  we take $`\mu ^{}=0.1`$. Selecting the values $`ϵ_M=1.8`$ and $`\lambda =0.3`$ and calculating $`T_c`$ with Eqs. (3-9), we obtain $`T_c=24.5`$K, which is close to the experimental value yamanaka . In absence of the plasmon part ($`\mathrm{\Gamma }_c^I=0`$) we obtain $`T_{c,ph}=0.5`$K which indeed confirms that the conventional electron-phonon mechanism cannot explain the high critical temperature observed in this material.

We point out that the calculation just performed for Li0.48(THF)yHfNCl makes use of reasonable, but still adjustable parameters $`\lambda `$ and $`ϵ_M`$. A more detailed analysis requires the experimental determination of these quantities prior to our calculation. It would thus be of interest to perform tunneling measurements which would allow to determine the function $`\alpha ^2(\mathrm{\Omega })F(\mathrm{\Omega })`$ \[$`F(\mathrm{\Omega })`$ is the phonon density of states whereas $`\alpha ^2(\mathrm{\Omega })`$ describes the coupling\], and correspondingly $`\lambda `$ (along with $`\mu ^{}`$; see, e.g., Refs. mcmillan ; grimvall ). Another method to determine $`\lambda `$ requires to measure the electronic heat capacity. Indeed, as is known, the Sommerfeld constant contains the renormalization factor $`1+\lambda `$ while the magnetic susceptibility is unrenormalized (see, e.g., Ref. grimvall ). Comparing these two quantities one can extract the value of the coupling constant $`\lambda `$. Such measurements, along with optical data would allow to carry out more detailed calculations of $`T_c`$ for specific metallochloronitrides.

In absence of such experimental data, we present in table 1 a few typical examples of calculated $`T_c`$ for various realistic values of the parameters $`\lambda `$ and $`ϵ_M`$ in Li0.48(THF)y HfNCl.

Note that in all cases the optical plasmon energy at $`𝐪=0`$ is of the order $`\omega _{pl,opt}(𝐪=0)11.3`$eV, in agreement with band structure calculation estimates. A more detailed analysis of other metallochlorinitrides will be described elsewhere.

In conclusion, the dynamical screening of the Coulomb interaction is an essential feature of layered structures that provides for an additional contribution to the pairing and leads to a drastic enhancement of $`T_c`$. The theory presented here enables us to give an explanation for the high critical temperatures observed in intercalated layered metallochloronitrides.