Monte-Carlo study of the Maier-Saupe model on square and triangle lattices

Physica 89A (1977) 522-538 © North-Holland Publishing Co.

MONTE-CARLO MODEL

STUDY

ON SQUARE

OF THE

MAIER-SAUPE

AND TRIANGLE

LATTICES*

RAYMOND D. MOUNTAIN? National Bureau of Standards. Washington, D.C. 20234 USA TH. W. RUIJGROK Instituut voor Theoretische Fysica, Rijksuniversiteit, Utrecht, The Netherlands~

Received 19 April 1977

The Monte-Carlo method is used to investigate the statistical physics of the Maier-Saupe model on square and triangle lattices. These systems are found to exhibit higher than first order phase transitions except for the negative coupling square lattice case. The transition is signaled by a pronounced maximum in the specific heat as a function of temperature. The mechanism of the transition is shown to be associated with the onset of partial ordering of the vectors in individual members of the ensemble. A temperature dependent internal order parameter is introduced to characterize the degree of internal order in the system since the ordering process does not break the rotational symmetry of the Hamiltonian and no macroscopic order parameter exists.

1. Introduction The Maier-Saupe

model

is e x t e n s i v e l y

u s e d as a m o d e l f o r liquid c r y s -

tals~'2). I n t h i s p a p e r w e d i s c u s s t h e r e s u l t s o f a M o n t e - C a r l o s i m u l a t i o n o f t h e e q u i l i b r i u m p r o p e r t i e s o f t h e m o d e l o n p l a n e l a t t i c e s . T h e m o d e l is d e f i n e d b y the Hamiltonian§

H =-JE

(t~i"/1i)2'

(1)

1,1

* Contribution of the National Bureau of Standards, not subject to copyright, partially supported by the Energy Research and Development Administration AT (49-16) 3003. + Guggenheim Fellow. Mailing address: Sorbonnelaan 4, Utrecht, The Netherlands. § The Maier-Saupe model is often expressed in terms of the Hamiltonian H = -J~

P2(r~i • ~j), i.)

where P2(x) is the second Legendre polynomial P2(x)= ~(3x: 1). This Hamiltonian is simply related to the one used in this paper by an energy shift and a rescaling of the coupling constant. The essential physics is unchanged by these shifts. 522

MONTE-CARLO STUDY OF THE MAIER-SAUPE MODEL

523

The sum is o v e r distinct nearest neighbor pairs on the lattice, J is the coupling constant, and t~i is a three dimensional unit vector fixed on site i. In the mean field approximation, the model exhibits a first order phase transition2). A Monte-Carlo study 3) of this model on a simple cubic lattice found a first order phase transition, albeit not at the t e m p e r a t u r e predicted by mean field theory. Our w o r k is concerned with the model on the square and triangular lattices and, in a g r e e m e n t with analytical studies4'5), we find no indication of a first order phase transition at non-zero temperature. We do, however, find evidence for a higher order transition which is associated with the d e v e l o p m e n t of "internal o r d e r " as the t e m p e r a t u r e is lowered. The "internal o r d e r " is not a macroscopically observable (thermodynamic) quantity, but the consequences of the order are readily observable in the specific heat. The plane lattice version of the M a i e r - S a u p e model is not expected to be appropriate for liquid crystals. It is an interesting model in its own right due to the opportunity it provides for studying a system with a higher than first order phase transition and for studying "internal order." At this point it is appropriate to indicate what is m e a n t by internal order. In the Monte-Carlo method (which is described in section 2), one constructs a sequence of configurations of the system which are felt to be representative of the thermal equilibrium condition. Properties of the system are constructed for selected configurations and the averages of these properties are equated with the canonical ensemble for the system. Internal order refers to order, present in individual m e m b e r s of the ensemble (configurations of the unit vectors t~i), which does not necessarily manifest itself in terms of a laboratory (external) f r a m e of reference. Internal order in the M a i e r - S a u p e model can be described in the following way. Consider the 3 × 3 matrix T ~° defined by 1

nTn~,

(2)

where n 7 is the a - c o m p o n e n t (a = x, y, z) of the unit vector at site i of the lattice. The sum runs o v e r the N sites of the lattice for a single configuration. Because of the overall rotational invariance of the Hamiltonian, the ensemble average of T ~° is e x p e c t e d to be 6"o/3. H o w e v e r , the eigenvalues of T ~° are not influenced by this invariance since they are internal to the s y s t e m for a given configuration. The ensemble average of the largest eigenvalue, A3, is the internal order p a r a m e t e r for the system. The orientation, relative to an external direction, of ~3, the eigenvector associated with the largest eigenvalue, is a measure of m a c r o s c o p i c order. For the s y s t e m s studied here, such m a c r o s c o p i c order does not exist. E v e n so, the internal order p a r a m e t e r (A3) has a well-defined t e m p e r a t u r e d e p e n d e n c e and provides a way of understanding the origin of the phase transition indicated by the a n o m a l y in the specific heat. A detailed discussion of the eigenvalues and eigenvectors of T ~ and of how they are related to more conventional concepts of order is found in the Appendix.

524

R.D.

MOUNTAIN

AND

TH.W.

RUIJGROK

The paper is organized as follows. Section 2 describes the Monte-Carlo simulation and the system properties constructed. The temperature dependence of equation of state and the internal order parameter are discussed for the square and triangular lattices for both the positive coupling (J > 0) and the negative coupling (J < 0) cases. The existence of higher than first order phase transitions is indicated by anomalies in the specific heat. Such anomalies occur in all but the square lattice-negative coupling case. Section 3 is devoted to the discussion of the origin of these anomalies in terms of the development of internal order.

2. Calculation The calculations were carried out using the Monte-Carlo simulation method introduced by Metropolis et al.6). Periodic boundary conditions were employed. The lattice arrangements studied are indicated in fig. 1 where 6 × 6 square and triangular lattices are illustrated. In the simulations, systems as large as 2 0 x 20 were studied. The dashed lines in fig. 1 enclose the actual system sites outside that boundary are taken to be periodic images of the enclosed sites. The different marks on the lattice sites denote the sublattice structures which are pertinent to the negative coupling (J < 0) studies. The Monte-Carlo algorithm introduced by Metropolis et al. 6) produces configurations appropriate to the Canonical ensemble. The quantities we have determined are; the mean energy per site, (E>, the fluctuations in the energy, which in the Canonical ensemble are related to the specific heat C by 7)

C/kB = ((/3AE)2),

(3)

the internal order parameter (A3) and the fluctuations in the order parameter ((AA3)2). In order to interpret the results of these computations, certain structural quantities were also constructed. The first of these is the probability +

+

+ O

+ 0

+

+

O

+

•

O

•

O

•

O

O

•

O

•

O

•

/-O

I+

+

0

0 0

O

+ 0

4-

0

O

I

O

•

O

•

O

•

O

+

I 'i 0

•

O

•

O

•

O

•

o

•

O

•

O

•

O

•

O

•

O

•

O

•

•

O

•

O

•

O

O

•

O

•

O

•

0

O

i

•

+ 0

O •

+

1. T h e

dashed

•

I I + 0 I I • + I - - .~ O

+ 0

•

+ O -- _ O__ "4-

Fig.

•

+

•

+ -- -- O

+

lines enclose

'4-

0

O

.2 O

q-

6 × 6 regions

of a square and a triangle lattice. The different

marks on the sites indicate sublanices important for negative coupling cases. The sites outside the d a s h e d l i n e s are t a k e n t o b e p e r i o d i c i m a g e s o f t h e e n c l o s e d

sites.

MONTE-CARLO STUDY OF THE MAIER-SAUPE MODEL

525

density of (r~i. 1~3)2 which m e a s u r e s the extent to which the principal eigenvector q~3 defines the orientation of the vectors on the lattice sites. Also constructed were densities of (ai" hi) 2 where (i, j) are either first or second neighbors on the lattice. These quantities provide a m e a s u r e of the degree of local order different f r o m that provided by the probability density of (hl • ~3)2. The probability densities, D ( c o s 2 O), were constructed as histograms using the rule indicated in eq. (4). The interval (0, 1] is divided into Im equal parts and if i = I/ Im, I = 1 . . . . . Ira, [ N u m b e r of o c c u r r e n c e s of i - 1 < cos 2 0 < i], of ] - a ~ s~-es]

O ( i ) A , = [ N u m b e r of ~ m - b e r

(4)

where the solid angle factor At is 27r A, = ~ [f)-

N/1 - I].

For our purposes adequate resolution was obtained with Im = 10, although Im = 25 was used occasionally. The calculations were carried out in the following way. A set of initial values for the x, y and z c o m p o n e n t s of the N unit vectors ni w e r e selected using a r a n d o m n u m b e r generator so that the vectors were distributed in a cone with a vertex angle of a b o u t 50 ° relative to the z axis. When more than one t h e r m o d y n a m i c state was simulated in one c o m p u t e r run, the final coordinates of run j - 1 were taken to be the initial coordinates of run j ( j > 2). Once initial coordinates and the t e m p e r a t u r e were assigned, the Monte-Carlo p r o c e d u r e was used to generate 200N trials in order to "equilib r a t e " the system. After the equilibration process was completed, 2 × 104N additional trials were conducted with a sample being taken e v e r y 2 N trials for a total of 10 4 samples per run. Tests were made which indicated that results are not influenced by the choice of initial conditions or by the n u m b e r of samples taken. The size of the system does influence the order p a r a m e t e r (h3), particularly in the transition region. The energy on the other hand is m u c h less sensitive to the size of the system. This is illustrated in fig. 2 where plots of (h3) vs. 1 I N are shown for several t e m p e r a t u r e s for the square lattice (J > 0). On this scale, plots of E vs. I [ N would be nearly horizontal straight lines. All calculations were p e r f o r m e d with J = _+~. This value was selected so that the /3 ~ ~ limit for the energy per lattice site for the square lattice with positive coupling would be ( E ) = - 1. This choice leads to ( E ) = _3 for the triangular lattice with positive coupling as fl ~ oo. The square lattice with positive coupling (J > 0) is the system we have studied most extensively and the description of it is the m o s t complete of those to be discussed in this paper. The qualitative features should be the same for nearly all plane lattices although some features do depend on the lattice structure. For this reason, the other systems considered were examined in much less detail.

526

R.D. M O U N T A I N

I

AND TH.W. RUIJGROK

I

I

0.9

6

,~ 0.8

.I- --j

~_ 07

~ o.0~ 0.5 0.4 I

i

0.0025

0 01

q04

1/N

Fig. 2. Size dependence of the order parameter f o r the square lattice (J > 0). Short calculations f o r a 40 x 40 lattice (fiN = 0.000625) suggest that the linear e x t r a p o l a t i o n f r o m I/N = 0.0l to I/N = 0.0025 g i v e s a c c u r a t e r e s u l t s for s m a l l e r v a l u e s of I/N.

2.1. Square lattice (J > O) The mean energy per site, (E), and the order parameter, (A3), are displayed as functions of/3 (inverse temperature) for the square lattice (J > 0) in figs. 3 and 4. These curves are based on the results for a 2 0 x 2 0 lattice. It is reasonable, considering fig. 2 and the associated discussion, to expect the corresponding curves for a larger system would show little change for the energy but would be somewhat changed for (A3). Most likely, (A3) would rise with increasing/3 more abruptly than is shown in fig. 4. It is noteworthy that the energy is a smooth function of /3 and gives no indication of a first order phase transition, in accord with theoretical arguments4). The low temperature 1.0

I

I

I

[

08

I

0.6

04

0

I

I

I

I

20

40

60

8.0

10 0

Fig. 3. M i n u s o n e t i m e s the m e a n e n e r g y p e r site for a 20 × 20 s q u a r e lattice ( J > 0) is s h o w n as a function of/3.

MONTE-CARLO STUDY OF THE MAIER-SAUPE MODEL 1.0

I

I

I

2.0

4.0

527

I

0.8

0.6

0.4 0

.

6.0

10,0

Fig. 4. The order parameter for a 20 x 20 square lattice (J > 0) is shown as a function of/3. b e h a v i o r o f t h e e n e r g y c a n b e d e s c r i b e d b y the e m p i r i c a l r e l a t i o n --(E} = 1.0-

1.4//3

(5)

f o r / 3 > 20. E q . (5) o v e r e s t i m a t e s - ( E } b y 0.025 a t / 3 = 10. T h e d i m e n s i o n l e s s s p e c i f i c h e a t C/ka is s h o w n in fig. 5 a n d t h e f l u c t u a t i o n s in t h e o r d e r p a r a m e t e r , ((AA3) 2) are s h o w n in fig. 6, b o t h for 20 x 20 l a t t i c e s . T h e s c a t t e r in t h e p o i n t s p r o v i d e s an i n d i c a t i o n o f t h e u n c e r t a i n t i e s i n v o l v e d I

I

I

I

I

O O

v

o

w

• °

..

4

I

I

I

•

5

6

7

8

B

Fig. 5. The specific heat for a 20 × 20 square lattice (J > 0) is shown as a function of /3 in the transition region. The scatter of the points indicates the uncertainty in the calculated values.

528

R.D. MOUNTAIN

T

l

--'r

AND

F

TH.W.

RUIJGROK

T

r

F

r

r -I

t

I 6

I

~ 8

I

4 ×

-A ,I.

v 2

0

°i 0

• •1 2

I

I 4

/3

10

Fig: 6. T h e m e a n square fluctuation in the order p a r a m e t e r for a 20 × 20 square lattice ( J > 0) is s h o w n a s a f u n c t i o n o f /3. T h e scatter in the points indicates the u n c e r t a i n t y in the calculated values.

in the determination of these quantities. S o m e states were studied using differing initial conditions. This had no effect on but changing the initial conditions frequently resulted in differing estimates for the fluctuations. This indicates that computer runs of much longer duration would be required to reduce the scatter in figs. 5 and 6. Such long runs were not undertaken as the available results are sufficient for our purposes.

I

I

I

I

I

O10t

005

\

o .4

B

8

X3 Fig. 7. D i s t r i b u t i o n s o f t h e v a l u e s o f A3 f o r a 20 × 20 square lattice ( J > 0) a r e s h o w n for three v a l u e s o f / 3 ; /3 = 4.1, 4.8 a n d 5.3: T h e a r r o w s m a r k the e n s e m b l e average Qt3). T h e s e distributions are t y p i c a l o f w h a t is f o u n d f o r states in the vicinity o f the transition region.

MONTE-CARLO STUDY OF THE MAIER-SAUPE

MODEL

529

The extent to which the square lattice is internally ordered undergoes a p r o f o u n d change as /3 increases f r o m 4.0 to 5.0. This is reflected in the fluctuations in (A3), which reach a m a x i m u m in the interval 4.6 < f3 < 4.9, and to a lesser extent the specific heat which also has a m a x i m u m in the same interval. These effects are a clear indication of a rapid but continuous variation in the internal order of the vectors on the lattice. Some insight into the type of changes in internal order which occur with increasing /3 can be obtained by examining the distribution of A3 values for which (A3) is the average and also to examine D[(ni" ~3)2], the distribution of the square of the cosine of the angle of individual vectors make with the principle eigenvector ~3. The distribution of A3 values for /3 < 4 . 4 is fairly narrow and peaked at (A3). As /3 increases beyond 4.4, the distribution broadens and the p e a k position b e c o m e s ill-defined. Starting with/3 = 5.0, the distribution is again narrow with the p e a k at (A3). Typical distributions are shown in fig. 7 for /3 = 4.1, /3 = 4.8 and /3 = 5.3. The arrow indicates the position of (A3). These distributions were constructed by sorting the A3-values with a resolution of 0.01 and normalizing the resulting distribution to unity. The line segments in fig. 7 are intended only as a guide to the eye. The density D[(a, • ~3)2] undergoes a dramatic change f o r / 3 on the order of 4.7-4.8. This is illustrated in fig. 8 where this density is shown for /3 = 4.5, /3 = 4.7, /3 = 4.75 and /3 = 4.8. The principal feature is the abrupt decrease in the fraction of vectors a p p r o x i m a t e l y normal to ~3 c o m p a r e d with the fraction roughly parallel to ~3. The absence of a preferred external direction was d e m o n s t r a t e d by constructing the distribution of values of the z - c o m p o n e n t of ~3 for states a b o v e and below the transition region. While the eigenvalues, A3, are fairly narrowly distributed (see fig. 7), the z - c o m p o n e n t s of 03 are widely distributed. At high

I

I

I

I

050

/

,_.>~ ' ~ 025

0.5

05

I

I

05

0.5

Fig. 8. H i s t o g r a m s showing the distribution of the orientation of vectors relative to the principal eigenvector $3 for four values of /3. T h e s e distributions for a square 20 × 20 lattice (J > 0) were c o n s t r u c t e d using eq. (4) with I m = 10. The distributions are normalized on a per-site basis.

530

R.D. MOUNTAIN AND TH.W. RUIJGROK I

I

i

40

3.0

"~=- 20

~=418I

/3 =4.1

1.0

/ / I

I

[1.5

0.5

/9=5.3

/ 0.5

Fig. 9. Histograms showing the distribution of the orientation of nearest neighbor pairs for a 20 x 20 square lattice (J > 0) for three values of/3. These distributions, which are normalized to 4-nearest neighbors, were constructed using eq. (4) with ,I,, = 10.

t e m p e r a t u r e s , the absolute values of z - c o m p o n e n t s o c c u r relatively f r e q u e n t l y in the 0 . 2 - 0 . 4 interval, i n d e p e n d e n t of initial conditions. At low t e m p e r a t u r e s , the z - c o m p o n e n t s are m o r e u n i f o r m l y distributed. T h e n e a r e s t n e i g h b o r and s e c o n d n e i g h b o r orientation correlation functions, the densities D[(fii • /~j)2] w h e r e (i, j) are nearest n e i g h b o r or s e c o n d n e i g h b o r pairs, p r o v i d e a n o t h e r w a y of e x a m i n i n g the internal o r d e r of the system. The c h a n g e s in these quantities are less d r a m a t i c than are the c h a n g e s o c c u r r i n g in D[(r~. ~3)2]. T h e first and s e c o n d n e i g h b o r orientation correlation functions f o r / 3 = 4.1, /3 = 4.8 and /3 = 5.3 are s h o w n in figs. 9 and 10, respectively. On the local level, the s y s t e m a p p e a r s to be s m o o t h l y and slowly e v o l v i n g t o w a r d s being o r d e r e d as /3 increases. T h e probability of finding first and s e c o n d n e i g h b o r pairs in a p p r o x i m a t e l y parallel (or antiparallel) configurations increases at the e x p e n s e of the energetically less f a v o r a b l e configurations w h i c h involve small values of (rJ/. ft.)2. T h e r e is no indication of the two p e a k e d structure f o u n d in D[(fi • q~3)2] f o r small values of/3.

MONTE-CARLO STUDY OF THE MAIER-SAUPE MODEL I

I

531

I

3.0

2.0--

~..=-

B =4.8 ].0

/ I

1

0.5

0.5

J

I 0.5

(6i .-~j)2 Fig. 10. Same as fig. 9 except for second neighbor pairs. 14

I

I

20

40

12

10

08

0.6

I

0

6.0

Fig. 11. Minus the mean energy per site for a triangle lattice (J > 0) is shown as a function of ft.

532

R.D. M O U N T A I N

I

0.8

AND TH.W. RUIJGROK

T

I

I

I

I

I

I

I

/~

-

m

A ,,<

v

06

i

n

04

m

I

20

4.0

I 6.0

Fig. 12. T h e o r d e r p a r a m e t e r for a 20 x 20 t r i a n g l e l a t t i c e (J > 0) is s h o w n as a f u n c t i o n of/3. N o t e n o u g h s t a t e s w e r e s t u d i e d to define a c u r v e ; q u a l i t a t i v e l y it is s i m i l a r to the s q u a r e lattice.

2.2. Triangular lattice (J > O) The energy per site, E and the order parameter, (A3), for the positive coupling case on a triangular lattice are shown in figs. 11 and 12. The results are qualitatively similar to those for the square lattice, the rapid change in internal ordering occurs o v e r the internal 2.9 < / 3 < 3.2. Another significant difference between the two lattice systems shows up in D[(~i. ~3)2], the probability density for the orientation of vectors relative to the principal eigenvector. This quantity for the triangular lattice does not exhibit the two p e a k structure at low values of /3 shown in fig. 8 for the square lattice. Instead, it is a smoothly increasing function of its argument. The roughly parallel orientations b e c o m e increasingly probable as /3 increases with the rate of increase being largest in the interval 3.05 < / 3 < 3.15. This is illustrated in fig. 13 for the region 2.75~
M O N T E - C A R L O S T U D Y OF T H E M A I E R - S A U P E M O D E L

533

0.75

/ 050

0.25

o

L 2.75

I 2.85

I 295

I 3.05

I 3.15

T 3.25

• 335

Fig. ]3. The distributions o f the orientation o f vectors relative to the principal eigenvector ~3 f o r the 20 x 20 triangle lattice (J > O) are shown f o r values of /3 which span the transition region. Each connected set of points corresponds to a specific value of I in eq. (4). The bottom line corresponds to I = 1 and the top line corresponds to I =Im = 10. The distributions are normalized on a per ion basis. The lines are intended only to guide the eye.

5.0

--

40--

30

--

J

<=20-

1.0 --

0

I

275

/ 2.85

I

I

T

•

295

3.05 B

3.15

325

335

Fig. 14. Same as fig. 13 for the distribution of the orientation of nearest neighbor pairs.

534

R.D. M O U N T A I N A N D T H . W . R U I J G R O K I 4.0

I

I

I

I

I

I

--

3.0--

" ~ - 2.0 c~

1.0--

0

I 2.75

I 2.85

I 295

I

1

~

I

3.05

3.15

3.25

335

Fig. 15. Same as fig. 14 for the distribution of second neighbor pairs.

2.3. Negative coupling The lowest energy state for positive coupling (J > 0) consists of all vectors aligned along a given direction. This state is highly degenerate* (2 u) because the Hamiltonian involves only (fii" fij)2. The lowest energy state for negative coupling (J < 0) is also highly degenerate but the details of the state depend on the lattice structure. The lowest energy is achieved when all nearest neighbor pairs are mutually perpendicular.

2.4. Square lattice (J < O) For the square lattice (see fig. 1) the lowest energy is achieved, for example, if the vectors on the " d a r k circle" lattice sites are oriented normal to the page and the vectors on the " o p e n circle" lattice sites lie in the plane of the page. Because the d e g e n e r a c y of this state is infinite, the tendency of the negative coupling square lattice to order both sublattices with decreasing t e m p e r a t u r e should be w e a k at best. This is bourne out by the energy and specific heat which are shown in fig. 16. At low temperatures, the specific heat reaches a plateau and nowhert; shows a peak.

* The degeneracy implied by the rotational invariance of the Hamiltonian applies to all states and is not considered specifically.

MONTE-CARLO STUDY OF THE MAIER-SAUPE MODEL

0.3

I

I

I

I 5

I 10

I 15

535

0.2 03 0

1.0 A ¢..4

v

0.5

o

20

Fig. 16. The mean energy per site and the specific heat for a 12x 12 square lattice (J < 0 ) are shown as functions of/3. The plateau in the specific heat for/3 > 10 is a reflection of the absence of kinematic constraints on the sublattices.

2.5. Triangular lattice (J < O) The lowest energy state for the negative coupling case (J < 0) on the triangular lattice is obtained when the vectors on the three sublattices indicated in fig. 1 are mutually perpendicular. This state has a slightly higher d e g e n e r a c y than the positive coupling case since the orientation of vectors on two sublattices can be interchanged without changing the energy. This extra d e g e n e r a c y is unlikely to seriously affect the properties of the system since the additional states in question are disjoint in phase space. The energy and specific heat for an 18 x 18 triangular lattice with negative coupling are shown in figs. 17 and 18. The energy changes rapidly o v e r the I

I

I

I

I

I

I

I

I

I

I

I

I

I

I

1

0.5

04

^ 03

0.2

0.1

0

I

2

I

4

I

/3

6

8

Fig, 17. The mean energy per site for a triangle Lattice (J < 0) is shown as a function of/3.

536

R.D. M O U N T A I N A N D TH.W. R U I J G R O K ,5.0

4.0

I

I

I

B

3.0 --

v

2.0--

10--

o

I

I

I

7.0

80

9.0

/3

Fig. 18. The specific heat for an 18 x 18 triangle lattice (J < 0) is s how n as a function of/ 3 in the transition region.

interval 7.75 < / 3 < 8.25 and this is reflected by the specific heat m a x i m u m at /3 - 8. The uncertainties in the energy fluctuation calculations preclude a more precise location of the m a x i m u m . In this same region, the principal eigenvectors for the three sublattices b e c o m e mutually perpendicular to a good approximation. The transition is a smooth one and the principal eigenvalues for each sublattice remain significantly less than unity even at considerably lower temperatures. For example, at/3 = 7.0, about 10% of the configurations generated have principal eigenvectors for the sublattices which are nonperpendicular by at least 20 °. At/3 = 8.25, less than 0.3% of the configurations were disordered to that extent. The principal eigenvalue for a sublattice is 0.46 at /3 = 7.0 and is 0.69 at /3 = 8.25. The fluctuations in the principal eigenvalue for a sublattice reach a m a x i m u m at about the same t e m p e r a t u r e as does the specific heat. In this respect, the positive and negative coupling cases show qualitatively similar behavior, although at differeot temperatures. The point to be stressed here is the interplay b e t w e e n the d e v e l o p m e n t of order within a sublattice and the d e v e l o p m e n t of order between the sublattices. The next section is devoted to the discussion of this interplay of ordering p r o c e s s e s and what it implies for the ordering in the positive coupling case.

MONTE-CARLO STUDY OF THE MAIER-SAUPE MODEL

537

3. Discussion

In this section we consider the way by which internal order develops with decreasing temperature in the plane lattice versions of the Maier-Saupe model. We begin by comparing the negative coupling cases for the square and triangular lattices. The triangular lattice version exhibits a peak in the specific heat while the square lattice version does not. This difference in the properties can be understood, as was mentioned above, in terms of the absence of strong kinematic constraints on the relative orientation of the sublattices in the square lattice model. Thus we are led to focus our attention on how internal order develops in the negative coupling, triangular lattice model. The most striking change which occurs in the transition region for the J < 0 triangular lattice model is the "lining up" of the principal eigenvectors of the three sublattices into mutually orthogonal arrangements. This is accompanied by increases in the principal eigenvalues of the individual sublattices which are similar to the change which occurs for the positive coupling case. The way in which the sublattices become ordered suggests the following picture of the events on the microscopic level for the positive coupling case. These events are reflected in the peak in the specific heat. At high temperatures the vectors are disordered. As the temperature is lowered some local ordering occurs without any significant spatial correlation of the locally ordered regions. The number of these regions increases with decreasing temperature and in the transition region, it becomes energetically favorable for these local regions to coalesce (the range of spatial correlation becomes large). This effect is marked by the peak in the specific heat which represents the onset of ordering of regions in the lattice. An examination of the coordinates of configurations for a large square lattice at/3 = 2.0, 4.0 and 6.0 confirms this picture of the transition. The behavior of the distributions for (ni " ~3)2 and for (r~i. t~j)2 for first and second neighbor pairs also supports this picture of the transition. The description of the ordering process proposed here is reminescent of Kadanoff's 8) way of introducing scaling into the study of critical phenomena in that above the transition region, the system consists of a number of small partially ordered regions whose size is scaled upwards during the transition. Of course, there are significant distinctions between our system and one which exhibits critical point behavior. The most important one is the absence of a first order phase transition at temperatures below the transition region; instead the coalesced, partially ordered ((A3)< 1) system becomes increasingly but continuously more ordered with decreasing temperature. Put another way, this system has a continuous s y m m e t r y (rotational invariance) which is not broken by the ordering process. We can state with certainty that the transition is not a first order phase transition. H o w e v e r , the limited accuracy of the computed values of the energy and the specific heat do not make it possible to determine if the transition is of second or higher order. This determination would require an analytic treatment of the transition.

538

R.D. MOUNTAIN AND TH.W. RUIJGROK

Acknowledgements T h i s w o r k w a s b e g u n w h i l e R D M w a s a g u e s t at the I n s t i t u u t v o o r T h e o r e t i s c h e F y s i c a , R i j k s u n i v e r s i t e i t , U t r e c h t . T h e w a r m h o s p i t a l i t y of the I n s t i t u t e is r e c a l l e d w i t h p l e a s u r e . T h e w o r k w a s s u p p o r t e d in p a r t b y the Division of Military Applications, Energy Research and Development A d m i n i s t r a t i o n , AT(49-16) 3003.

Appendix T h e i n t e r n a l o r d e r p a r a m e t e r , (A3), is t h e e n s e m b l e a v e r a g e o f the l a r g e s t e i g e n v a l u e o f t h e m a t r i x T ~ i n t r o d u c e d in eq. (2). T ~ is a real, s y m m e t r i c , p o s i t i v e definite, 3 x 3 m a t r i x w i t h d i a g o n a l e l e m e n t s w h i c h lie in the i n t e r v a l 0
and

T~/°

= hi~

( a , / 3 = x, y, z).

T h e s e p r o p e r t i e s a r e e a s i l y p r o v e d b y o b s e r v i n g t h a t T ~ = 1 a n d f ~ T ~ t ~ f t3 >1 0 f o r a n y f. In a n a l o g y w i t h t h e rigid b o d y problem9), the e i g e n v a l u e s c o r r e s p o n d to t h e p r i n c i p a l m o m e n t s o f i n e r t i a o f t h e c o l l e c t i o n o f v e c t o r s on t h e l a t t i c e a n d the e i g e n v e c t o r s d e t e r m i n e t h e o r i e n t a t i o n of t h e i n e r t i a e l l i p s o i d . In this w a y , A3 is a m e a s u r e o f the e x t e n t to w h i c h i n d i v i d u a l c o n f i g u r a t i o n s d e p a r t f r o m an i s o t r o p i c one. A s s u c h (A.3) is an i n t e r n a l o r d e r p a r a m e t e r . T h e u s u a l o r d e r p a r a m e t e r e m p l o y e d in the d e s c r i p t i o n o f liquid c r y s t a l s 2) w o u l d be o b t a i n e d b y first a v e r a g i n g T ~° a n d t h e n finding t h e e i g e n v a l u e s . F o r this s y s t e m w h i c h p r e s e r v e s t h e r o t a t i o n a l i n v a r i a n c e o f the H a m i l t o n i a n , the u n i n t e r e s t i n g r e s u l t hi = A2 = A 3 - - ~ w o u l d r e s u l t f r o m this o p e r a t i o n .

References 1) P.G. deGennes, The Physics of Liquid Crystals, (Oxford Univ. Press, London, 1974). 2) M.J. Stephen and J.P. Straley, Rev. Mod. Phys. 46 (1974) 617. This review article is a good source of references. 3) P.A. Lebwohl and G. Lasher, Phys. Rev. A 6 (1972) 426. 4) P.A. Vuillermot and M.V. Romerio, Commun. Math. Phys. 41 (1975) 281. 5) H. Kimura, J. Phys. Soc. Japan 37 (1974) 1204. 6) N. Metropolis, A.W. Rosenbluth, M.N. Rosenbluth, A.H. Teller and E. Teller, J. Chem. Phys. 21 (1953) 1087. 7) T.L. Hill, Statistical Mechanics, Principles and Selected Applications (McGraw-Hill, New York, 1956)p. 101. 8) L.P. Kadanoff, Physics 2 (1966) 263. 9) H. Goldstein, Classical Mechanics (Addison Wesley Publ. Co., Reading, Mass. 1950) Ch. 5.