Thermodynamics of an elastic ferroelectric model

Thermodynamics of an elastic ferroelectric model

Physica 67 (1973) 291-307 0 North-Holland PubIishing Co. THERMODYNAMICS FERROELECTRIC A. HINTERMANN OF AN ELASTIC MODEL and S. P. OHANESSIAN Lubo...

839KB Sizes 1 Downloads 40 Views

Physica 67 (1973) 291-307

0 North-Holland PubIishing Co.

THERMODYNAMICS FERROELECTRIC A. HINTERMANN

OF AN ELASTIC MODEL

and S. P. OHANESSIAN

Luboratoire de Physique Thiorique, Ecole Polytechnique FJdtGale, Lausanne, Suisse

Received 7 August 1972

synopsis The modified KDP model of Wu on a compressible harmonic lattice is considered. An electric field and normal forces act on the system. The specific heats, mean lattice spacing, polarization and susceptibility are computed. Elasticity does not change the order of the transition, which remains of second order (without latent heat). The specific heat at constant volume is finite through the critical temperature, while the speciiic heat at constant force diverges at T,+. The dependence of the critical temperature on electric and force fields as well as the isotherms are discussed.

1. Introduction. Numerous models defined on rigid lattices have been proposed for a description of systems undergoing a phase transition, for example the Ising model for the order-disorder transition in binary alloys and for the magnetic order in strongly anisotropic insulators. Other models such as the potassium dihydrogen phosphate (KDP) model of Slater and the F model of Rys describe ferro- and antiferroelectrics. However, all these models do not take into account the displacements of the atoms from their equilibrium positions. Several authorslW4) have introduced elasticity in these models, by considering explicitly a dependence of the interactions on the mean lattice spacing. More realistically, other authors5-lo) attributed this dependence to the local lattice spacing. The main idea is to reduce the elastic problem to that of a rigid one. First to introduce elasticity in a ferroelectric KDP model was Coplang). His work, published in letter form, starts with a hamiltonian in terms of displacements and gives the partition function and specific heat at constant volume. The present paper is based essentially on his work, giving details of the calculations, and studying more extensively the thermodynamic properties of the system, especially the electric ones. The difference between Coplan’s work and ours is in the choice of the electric energy. This results in an important difference in the behaviour of polarizability and electric susceptibility of the two models12). 291

A. HINTERMANN

292

AND S. P. OHANESSIAN

In section 2, we start from the hamiltonian of the rigid model, and we introduce elasticity following the procedure used in other articles5*6*8). The lattice will be treated in the harmonic approximation. In section 3, we reduce the partition function of the elastic model to the partition function of a rigid model in the presence of an effective field. The differences between our partition function and the one of Coplan is due to our modification of the electric energy and some modified definitions of his constants. The thermodynamics of the system in the frame of the modified rigid KDP model of Wu”) is studied in section 4. We first recall the known results of the modified KDP model, then compute the mean lattice spacing, the specific heats and the electric properties. Since in this rigid KDP model the sign of the polarization cannot be changed, it describes a pyroelectric rather than a ferroelectric system. We study then the dependence of the critical temperature on the external fields and deduce constraints for the elastic constants. Furthermore, we investigate the equation of state, and conclude that harmonicity does not change the order of the transition, which remains of second order. Finally we show that the elastic model , exhibits a piezoelectric behaviour. 2. Hamiltonian and phase space. The system is defined on a lattice generated by two fundamental translation vectors a, and a,,. There are N = (2n + 1)’ lattice points, called vertices, and labelled by (i, j) where -n < i, j < n. A primitive cell comprises a cation labelled by (i, j) and two anions labelled by (i + 3, j) and (i, j + 4).

0

:

occupied

Fig. 1. Equilibrium

positions

anion position

0

: vacant anion position

of the anions surrounding the cation d 1.3 = 4, = -1, uy,j = u+J = 1).

(i,i) (configuration

with

THERMODYNAMICS

OF AN ELASTIC

FERROELECTRIC

MODEL

293

We call “rigid system” the model where the cations (i, j) are at vertices (i, j) and the anions can be found on one of the two equilibrium positions on the lattice bonds (fig. 1). We thus have 24 = 16 possible configurations for the anion positions surrounding the vertex (i, j). In order to characterize a rigid configuration, we associate with each cation (i,j) four spin-like variables u:,~ (y = r, 1, u, d for right, left, up, down). If the (i,j, y) anion is located at the equilibrium position at (i,j) E aY_ (resp. a?+), we define a:,, = -1 (resp. + 1). The a:,, and u: are not all independent. Obviously we have the identities: Uf.l = -a:,,+1 a’_ + a:

and

oY**= -c;+1.*,

= a; + a!_ = I&..,

a”_ + ad, = a:

-n
(1) (2)

Rigid configurations could also be defined by means of an arrow pointing toward (resp. away from) the vertex (i, j) if u:,~ = - 1 (resp. + l), see fig. 3. For each of the 16 arrow configurations of a vertex, there is an energy assignment .Q ((a;,,}) and the problem is to compute the partition function:

We define now an “elastic model” where all ions can move. The system is in interaction with a normal force and an external constant electric field. The energy asstgnments et, I are now functions of the displacements {Ri, I} of the anions relative to the cations. At a fixed vertex only 6 of the 8 components of R:,,are independent. Thus the hamiltonian reads: H = 1 H,,, + work of force + electric energy, 1,1 H,,j = K-E. + ei,j ({~Y.JPRi,,}),

(4)

where K.E. is the kinetic energy of the cations and anions. We shall now discuss the terms of this hamiltonian. We expand cl, J around 0:; I, the equilibrium positions of the system in the ’ presence of the fields:

where 01and /I are Cartesian coordinates. Here, we assumed that el, j depends only on the displacements of the four neighbouring anions relative to the cation (i, j).

294

A. HINTERMANN

AND S. P. OHANESSIAN

This implies having small Boltzmann factors for large displacements compared to the lattice spacing, i.e., requiring large A$’ compared to E~,~({u~,~))/u~.In this way we can associate with each elastic configuration a unique rigid one. The two first terms of (5) are the harmonic terms of the expansion. In order to have an exact and easily solvable model, we choose appropriate A$ such that E~,~reads: &i,j

({Ur,j, Rr,j}) =

Ei,

j

{Kr (ET,,/- ar)z + sy(4~~.i)2~

({Uzj}) + C + Y

1 +o:, ( +c Y 2

K: (ET,J - a:)’ )

where 5 and 7 are defined by: R:,, = (51..,, r:.,),

Rt,j = -(6:.~, qf.,),

Rll, = (r:,j> E:,,),

&

= -($,r,

&).

In order to determine a point of the phase space of the elastic model, we have six momenta, a priori eight displacement coordinates, and the corresponding equilibrium configuration. Thus, provided the hypothesis implied by eq. (5) is correct, we have for the phase space r, I’ = R6N @I R8N @IcrBzN, where R6N, RaN, 0@2N = {-1,1}@‘2N are respectively, the momentum, displacement, and spin subspaces. In fact, there are only six displacement degrees of freedom per vertex, and so appropriate projections of the RsN space on R6N subspaces are performed.

Fig. 2. An elastic configuration

characterized

by Ri,i.

THERMODYNAMICS

OF AN ELASTlC

FERROELECTRIC

Finally, the work of the forces W, is proportional of a cell, and the electric energy W, is proportional fore

MODEL

295

to the variation, of the length to the dipole moment. There-

The total hamiltonian with ai,, ({&, R:,,}) written as in eq. (6) is the starting point of the work of Coplang). The difference between Coplan’s hamiltonian and (4) consists in the form of the electric energy. Coplan considers terms up to quadratic order in the displacements, which leads to a polarization tending to infinity for high temperature12). Moreover, his linear terms do not contain the total dipole moment; indeed, choosing appropriately his constants, one has:

which is different from E * Rr, k. 3. The partition function of the elastic model. There are N = (2n + 1)2 vertices, and 12 degrees of freedom per bulk vertex. Taking into account the surface, we have a total of 12 (2n + 1)’ - 8 (212+ 1) = 4~ degrees of freedom. Thus the partition function reads :

QN = g

_;Idi _rCd”0(d’d ew ( -DO

= 2 &({a>),

(8)

where m is of the order of the lattice spacing. Within the hypothesis contained in (5), the domain of integration can be extended to infinity. The integration over the p’s is immediate, thus giving: (A,)-‘” (&)-8n(2n+ l), where A,,,, = (p/2xmv, h)3. We define:

For a given spin configuration all the variables in the hamiltonian dent. Thus for ui,, = + 1 ( - 1) we have the contribution :

are indepen-

(a:,, = +I):

(n/B>* G: =

rCd&) exp { -B Wi (5:~ - a;)’ + t~.,Fxll,

-co

G; = (K$)-* exp /I [(F,/4K;)

- Fxa;],

296

A. HINTERMANN

AND S.P. OHANESSIAN

/3 (I;x + Q2

G’_ = (K:S)-‘exp

4Kf.

E,” + 4s’ - (Fx + Ex) a’_ . >

The other G$ are defined by the convention: E-t-E;

Fx + Fx',,

E, t-) E,;

F, + F,,

E,t,

-E,;

Thus for a spin configuration

F,+F,.

{o> we have:

&CM>= ,2(2n+1)31;2~~~8~(2~+l)(X/~)12nZ+6n

The factor m2(*“+l) can be disregarded in the thermodynamic tion function will then be:

QN = m

2(2n+1)~;2N3iban(2n+l)

x

limit. The parti-

(,g)12n2+6n

(13)

where

and

is the partition function of a rigid model in the presence of the effective field w. Having adopted the same notation as Coplang), his partition function and ours are formally the same. The differences in the expressions of the G” are due to our different choice of the expression of the electric energy. Even for E, = E,, where the linear part of Coplan’s electric energy (7a) coincides with ours, the

THERMODYNAMICS

OF AN ELASTIC

FERROELECTRIC

297

MODEL

electric properties of the system cannot be identical, since:

lim wGf(E,,

E,. E,-‘E

E,J f

>

Af(E,E), -&f(E,

E) .

Let us note that the result is valid without restrictions on the {c} [e.g., as Coplan pointed out, the ice condition is not necessary to deduce eq. (13)J. 4. Thermodynamics of the model. The rigid model we have is generally called a “16 vertex model”13). The study of the thermodynamics of the elastic model requires a solvable rigid model in the presence of a field [see eq. (13)]; the only model which has been solved in presence of an electric field is the modified KDP model of Wu”). Like Coplan, we restrict the analysis of the thermodynamics to this model. First we shall briefly recall the main properties of the model. 4.1. The modified KDP model. The modified KDP model issues from the general 6-vertex Slater KDP model with two main restrictions: eI = co, e3 + e, = e5 + e6. Thus, the modified KDP model is in fact a 5-vertex model, since the first configuration is forbidden. Wu obtained the solution by using the pfaffian

SLATEI? wli

Do

Fig. 3. Vertex-energy

I,

assignments

8,

I

I

‘I

for the Slater KDP model and the modified one in presence of an electric field 8.

method, and gave the expressions of the free energy, the internal energy and the polarization. The condition for criticality depends on 8’ and divides the 8 plane into three regions (fig. 4). This is due to the fact that the ground state is different in each region. On the boundaries, the critical temperature T, is zero. This model has a second-order phase transition (no latent heat). The specific heat C diverges as (T - T,)-’ for T > T,. Below T, the configurational restrictions (ice condition) freeze the system in the configuration e2 (e 3, e4) of fig. 3, so the polarization of the model is equal to the ground-state polarization and C = 0. The freezing of the system is independent of the energy assignments E and the electric field, and will therefore persist for the elastic model.

298

A. HINTERMANN

AND S.P. OHANESSIAN

Fig. 4. Three regions in the 8 plane and the corresponding

ground-state

configuration.

4.2. Thermodynamic functions. Since the configurational energies of the four last vertices (e3 = e4 = e5 = e 6 = E) in the Wu model are equal, we impose that the elastic energies of these configurations are still equal. This implies a reduction of the number of elastic constants. The only possibility is to identify

(K”+,4, S”) =

(K; , a;,

9)

and

(K:, ai, S’) = (Kg,ad,,

Sd),

which imply laxI = J(L,!. From eq. (13), we obtain the expression of the Gibbs function 9 in the thermodynamic limit by the relation:

We now define the following quantities :

THERMODYNAMICS

OF AN ELASTIC

FERROELECTRIC

MODEL

299

Thus we have #%P = 2 In I, + 4 In jlh - 3 In (7r/b) + $ In (K:Kf.S’K:K!_S’)

- 38 [tF2x+ + &E2 (cr+ + t+)

- 21;. a - +E. Fu_ + (Ex + E,) d-1 + ,W”“’ (,!3,y) = /W”“(/?,

E,F) + jWrig(~,~)>,

where Rrtg is the free energy per vertex of the Wu rigid model and a = (a, a). In order to get thermodynamic quantities at constant volume, we make a Legendre transformation, taking the canonical conjugate of F. From the hamiltonian, it follows that the canonical conjugate of FXFx(,,) is the mean distance IXo,)between two nearest vertices (iW/aF),,. = 1.We can now define a new Gibbs function ~(B,E,~)by: g(B,E,Q

=~((B,E,F(B,KO)

- F(B,E,O*I.

(16)

In the following, we restrict ourselves to E, = E,, = E and Fx = F,, = F implying yX = yY = y. These points are located either in region I or on the boundary of II and III (fig. 4) where T, is zero. To have a nonzero T,, we have to restrict C, E, F in such a way that ?,u< 4s (see 4.5). 4.3. Thermal quantities. In this section, we calculate the mean lattice spacing I, = I, and the specific heats at constant force C, and at constant volume CI. For the mean spacing we have:

zi =

(g>,,, = (z),,, +(F),(-$),,:

where Prig= -cWrig/i3y is the polarization prig =

(144

of the rigid model,

((2/n) arc cos U - 1, (2/n) arc cos U - 1) [ (-1,

i =x,Y T > T, T<

-1)

T,

and U = Cexp ( -@ [&F’x_ + $E2 (CJ- + t_) - +EFa+ + Ed+ - E]}. The critical temperature

Tcll) is determined by U = 1. We thus obtain:

1, = a - $F+C+ + fEi~_ - *PIi’ (Fix- - EiO+). For the computation

(154

(If%

of the specific heats, we shall work with E = 0 and F,

= Fy = F.

6 + 2p2U(&_F2 C, -= k,

- E)~

x (1 - .2)+ 6

T>

T,, (17)

T<

T,.

300

A. HINTERMANN

AND S. P. OHANESSIAN

Fig. 5. Specific heats versus inverse temperature.

Thus, the specific heat at constant force diverges as (T - T,)-’ as in the rigid model of Wu. In order to evaluate the specific heat at constant volume C,, we calculate first the energy % (b, 2) from the 9 (/?, I) function.

=

aw (AagF (A 0)>,

+

(““yy))p(g), - (Y), (18)

-1-F.

The specific heat Ci is then:

Cl -= k,

I

Lx+ - x_ (1 - 2x-l arc cos U)] (*F2x_ - a)2 x Lx+ -x_(l

-2x-‘arccosU)](l

- U2)f+x-1@2F2~:U

(19)

T > T,

I -

6

T < T,

Since we have chosen Fx = F,,, which implies in particular that I, = 1, = 1,and since we have forbidden tangential constraints, a derivative at constant volume is

THERMODYNAMICS

OF AN ELASTIC

FERROELECTRIC

301

MODEL

identical to a derivative at constant 1. We note that for the elastic model C1 does not diverge (except for F = 0), while C, and the rigid-model specific heat do. 4.4. Electric properties. We shall consider here the polarization, easily get the electric susceptibility. We have by definition

and thus

P = -ac9 (b, E, F)/aEIF

and we easily obtain 4P, = Ei (u+ + t+) - F~o_ + 26_ + Pf” (EGO- - Fig+ + 26,) + P;“E~E_.

01

0

01

0.2

03

04

Fig. 6. Polarization

/$c

q

(20)

0.64

versus inverse temperature.

For the susceptibility xu = aPi/aEj we have c+ + t+ + Pp (6_ + t-)

4Xli =

I

DU + (1 -

L(o++ t+) B

-T-

4XYX= 4&

=

x

uq+

((Eia-

- (a_ + t_)

- Fia+ + 26,)’ T<

+ EFt?)

(21)

T,

u [2E,E,o_ (1 - u*)+ - (ExF, + E,Fx) c+ + 2 (Ex + Ey)

0

T >To

T<

T,.

d+l

T> Tc,

302

A. HINTERMANN

AND S. P. OHANESSIAN

We see that the polarization is continuous at the critical temperature, susceptibility exhibits a divergence at Tz. 4.5. The critical

temperature.

while the

The critical temperature given by

U(T,*, F*, E*) = 1 is determined by C + F*2 + alEe

-T,*ln

- ol,E”F* + E” -

1 = 0,

(22)

where we have used the dimensionless quantities T* = (/!I&)-‘,

F* = 4 (X-/E)+ F,

E* = (6+/a) E,

(23) (x1 = (a- + r-) E/48:,

a2

=

(a+P+>

The parameter 6 = 01~- && determines the With 6 positive (resp. zero, negative) eq. (22) hyperbolic) paraboloid (fig. 7). As pointed out in section 4.2, the effective smaller than 3 in order to have a nonzero

(4d

shape of the critical surface U = 1. describes an elliptic (resp. parabolic, field y* (T*, E*, F*) E Y/E must be critical temperature. Therefore we

T'

1 ’

\

1’

\

r’,_----‘?_ / . ‘.

\

‘.

_

I

E'

+ a)

6 =

0.0075

c

0.37

=

Fig. 7. The shapes of the critical temperature surface. The dash-dotted line indicates the line of the maxima of the parabola at E* = const. and their projection on the Z’* = 0 plane.

should also study the surface p* = 3 which, in terms of the dimensionless quantities introduced by eq. (23), is described by -T*

In 2C + Fe2 + ,x,E*~ - a2E*F” + E* = 1.

(24)

THERMODYNAMICS

OF AN ELASTIC

FERROELECTRIC

MODEL

303k

The two surfaces U = 1 and 2y,* = 1 intersect in the plane T* = 0 and the first is obtained from the latter by a scaling of In C/in 2C of the T* axis. Let us look for the intersections of the two surfaces with a plane E* = const. The shape of such intersections is independant of the value of 6. In fig. 8 we illustrate the possible cases.

/-A

Y’- %

a)

O
‘“\// I/ -

!4=’r

\A

b) 1/2-=c < 1

c)

C>l

Fig. 8. Intersections of the surfaces 2~* = 1 and U = 1 with planes E* = const. and the variation of the effective electric field y* as functions of temperature in the d place (see fig. 4).

a) 0 < C < 3. For F* = F? the effective field at zero temperature is located in region I of fig. 4. If we increase the temperature, y*(T,*) is still in region I and we always have y* cc T*. For F* = F,*, y* (T* = 0) lies on the intersection of regions II and III. Therefore the ground state of the system is given by any allowed arrangement of configurations e3, . . . , e6 and the system will be in a paraelectric phase for all temperatures. b) 3 < C < 1. For F* = F? the effective field is located for all temperatures in region I. The system undergoes a phase transition and we have y* cc -T*. For F* = Fz, y* (T* = 0) starts in the paraelectric phase and enters with increasing temperature in region I without being critical for any temperature. c) C > 1. There exists no critical temperature for F* = FT and the system is for all temperatures in a completely polarized state. If F* = F! the system is in a paraelectric phase at T* = 0, undergoes a “phase transition” at T,* and will be in a completely polarized state for high temperatures.

A. HINTERMANN AND S. P. OHANESSIAN

304

Now we shall restrict ourselves to the case C < 1, the other values of C leading to noninteresting models. 4.6. The equation of state. Let us consider a system with given elastic constants. At fixed E*,F* and T* the equation of state, eq. (16 a), depends on whether T* is below or above the critical temperature. If T* < T,*the “spin” part of the model is in a completely polarized state with Prig = ( - 1, - l), and the equation of state reads, with I* = l/a and 01~= (~+/2a) (E/X__)‘, I* = 1 + &m2cx6E"(y- 1) + cw,F*(y- l),

T* < T,*,

(25)

where we have set y = S-/S, = a-/o+ = x-/x+ = z-/r+ in order to reduce the number of elastic constants. Thus below the critical temperature the system is governed by the harmonic part of the hamiltonian and therefore the mean lattice spacing does not depend on the temperature. For temperatures T* > T,*, the thermodynamics of the system is dominated by the “spin” part of the hamiltonian. The mean lattice spacing depends on the temperature uiu the rigid polarization Prig and reads 1” = 1 + j-p2a6E*(y + Prig) - R~F* (1 + yPri8), T* > T,*,

(26)

where Prig= (2/x) arc cos U - 1 and U = Cexp[-b*

(F*2 + n,E*2 - (x,E*F*+ E* - l)].

The surface of critical temperatures in the F*, E*, T* space separates the domains of validity of the two equations of state. Now we shall discuss the isotherms in the I*, F*,E* space. Let us draw the isotherm at temperature Tf in the case 6 > 0 and C < 3. If T: < T,*(F*,E*) the isotherm T: contains an elliptic part of the plane described by eq. (25). Since (8Z*/dF*)E.,Ts < 0 and (8Z*/dF*)E*,T,e = -00, the isotherm takes the form indicated in fig. 9. We have excluded models with C > 1 since they are not physical. But even the cases C < 1 are not without ambiguities. i) All isotherms collapse in a part of the plane of eq. (25). This is due to the fact that the constant term of eq. (25) does not depend on the temperature. ii) Due to the rather special temperature dependence of the equation of state, eq. (26), all isotherms pass through a line. Indeed, if F* = +x,E*, eq. (26) is independent of Prigand thus its temperature dependence drops out. Moreover, if

THERMODYNAMICS

OF AN ELASTIC FERROELECTRIC

MODEL

305

Fig. 9. The isotherm T;’ = 20 in the case 6 = 0.0075. The dotted lines represent I* on planes of constant E*. The interior of the smaller ellipse is the domain of validity of eq. (25) (T,* s T$).

F*= + cx,E*= - oc,E*F* + E* -

1 = 0, we have U = C; therefore I* is again independent of the temperature, thus giving a second-order curve in the I*, E*, F* space where all isotherms intersect (fig. 9). Besides these curves, different isotherms do not intersect, and behave as illustrated in fig. 10. Wu has pointed out that the phase transition in his modified KDP model is of second order (without latent heat); while in the Slater KDP model, Lieb14) has shown that for zero field the transition is of first order (with latent heat). However, the specific heat of the two models exhibits the same behaviour through T,.

Fig. 10. The isotherms

T$ < T$ < T: for E* = const.

It follows from the shape of the isotherms (see fig. 9) thqt our elastic model always undergoes a second-order phase transition since the curve of critical temperatures coincides with the limit of stability. Thus the elastic part of the hamiltonian does not change drastically the rigid-model behaviour. WagneF) obtained a similar result for the harmonic Ising model.

306

A. HINTERMANN

AND S. P. OHANESSIAN

5. Conclusions. Elasticity has been introduced in the KDP model in a phenomenological way; i.e., by expanding the energy assignments in the relative displacements. The fact that we restricted the expansion to harmonic terms has important implications which we want to discuss. i) Obviously the system is unstable with respect to tangential constraints. ii) By an appropriate choice of the fields F and E, the system can collapse, namely the mean lattice spacing can be zero. iii) Different isotherms have common parts in the E*, F*, I* space. iv) Elasticity introduced in a harmonic way does not change the nature of the transition, which remains that of the rigid model (second order, without latent heat, for the modified KDP model of Wu; first order for the Slater KDP model). v) Above the critical temperature, the system is governed by the rigid model, whereas below T,, the harmonic part dominates. The only effect due to the harmonic part is a smearing out of the singularities at T, for thermodynamic quantities at constant volume (e.g., C,). The instability with respect to tangential constraints is suppressed if we consider anharmonic terms in the expansion of the energy assignments. Moreover, they lead to a = a(T), and, therefore, different isotherms have no common part; thus the ambiguity of intersecting isotherms will be lifted. As in the elastic Ising-model case 5,8*10)we expect that anharmonicity leads to a first-order phase transition. The model exhibits a piezoelectric behaviour. Increasing the force at constant temperature as indicated by the pathr of fig. 8, we note that the system goes from a completely polarized state to a partially polarized one. The same argument holds when we vary the electric field. Finally, the system is by definition anisotropic since the first configuration e, (see fig. 3) is forbidden. Therefore an isotropic limit (all elastic constants equal) makes no sense. Acknowledgments. We are indebted to Professors Ph. Choquard, A. Quattropani, G. H. Wannier, and particularly Dr. Ch. Gruber, for valuable discussions and suggestions. Thanks are due to Mlle. M. Pahud for programming work.

REFERENCES 1) 2) 3) 4) 5) 6) 7) 8) 9) 10)

Rice, O.K., J. them. Phys. 22 (1953) 1535. Domb, C., J. them. Phys. 25 (1956) 783. Mattis, D.C. and Schultz, T.D., Phys. Rev. 129 (1963) 175. Garland, C.V. and Renard, R., J. them. Phys. 44 (1965) 1130. Baker, G.A. and Essam, J.W., Phys. Rev. Letters 24 (1970) 447. Wagner, H. and Swift, J., Z. Phys. 239 (1970) 182. Coplan, L.A. and Dresden, M., Phys. Rev. Letters 25 (1970) 785. Jasnow, D. and Wagner, H., Z. Phys. 249 (1971) 101. Coplan, L.A., Phys. Letters 35A (1971) 309. Horner, H., Z. Phys. 251 (1972) 202.

THERMODYNAMICS

OF AN ELASTIC

FERROELECTRIC

MODEL

307

11) Wu, F. Y., Phys. Rev. Letters 18 (1967) 605. 12) Favre, P.P., diploma thesis EPF-L (1971), unpublished. 13) Lieb, E.H., Lectures in Theor. Phys. XID (1968) 329, K.H.Mahanthappa and W.E.Britten, eds., Gordon and Breach (New York, 1969). Lieb, E.H. and Wu, F. Y., Two Dimensional Ferroelectric Models in Phase Transitions and Critical Phenomena, C.Domb and M.S.Green, eds., Acad. Press (New York, London, 1972). 14) Lieb, E.H., Phys. Rev. Letters 18 (1967) 108.