Thermodynamic properties of ferroelectric NH3CH2COOH·H2PO3 crystal

Physica B 520 (2017) 164–173

Contents lists available at ScienceDirect

Physica B journal homepage: www.elsevier.com/locate/physb

Thermodynamic properties of ferroelectric NH3CH2COOH·H2PO3 crystal a

b,⁎

b

I.R. Zachek , Ya. Shchur , R.R. Levitskii , A.S. Vdovych a b

MARK

b

Lviv Polytechnic National University, 12 Bandery str., 79013 L'viv, Ukraine Institute for Condensed Matter Physics, 1 Svientsitskii str., 79011 L'viv, Ukraine

A R T I C L E I N F O

A BS T RAC T

PACS: 77.84.-s 64.60.Cn 77.22.-d 77.80.-e 77.80.Bh 77.65.Bn

Using a modified microscopic model of NH3CH2COOH·H2PO3 by taking into account piezoelectric coupling with strains εi, ε4, ε5 and ε6 in two-particle cluster approximation, the temperature dependence of polarization and tensor of static dielectric permittivity of mechanically clamped and free crystal, their piezoelectric characteristics, elastic constants and heat capacity are calculated.

Keywords: Ferroelectricity Spontaneous polarization Dielectric permittivity Heat capacity Hydrogen bond

1. Introduction Glycinium hydrogenphosphite crystal, NH3CH2COOH·H2PO3 (GPI), is a very interesting compound due to the combination of structural elements typical of different classes of crystals. On the one hand, it contains the covalently bonded phosphite HPO3 groups which are linked through the hydrogen O − H ···O bonds thus forming the chains running approximately along the c-axis (see Fig. 1). Such structural components are usual for non-organic ferroelectric materials, in particular for crystals of KH2PO4 family. Exactly the hydrogen bonded PO4 groups play a leading role in the mechanism of ferroelectric phase transition in KH2PO4-type crystals. On the other hand, there are four organic glycinium groups NH3CH2COOH in GPI unit cell which are linked by four additional hydrogen bonds with phosphite HPO3 groups belonging to two different phosphite chains (see Fig. 1). Such a structural complexity may result in great diversity of physical properties of this material. Therefore, GPI crystal has been intensively studied by means of dielectric [1,2], acoustic [3–5], calorimetric [6], optic [7] and spectroscopic [8,9] methods. A theoretical description of dielectric properties of GPI was presented in [10,11]. At room temperature, GPI crystalizes in a monoclinic P2 1/ a space group (No. 14, Z = 4) which transforms to P21 symmetry (No. 4, Z = 4) [12–14] below structural phase transition occurring near TC = 224 K. High-temperature phase is non-polar paraelectric phase whereas the

⁎

Corresponding author. E-mail address: [email protected] (Y. Shchur).

http://dx.doi.org/10.1016/j.physb.2017.06.013 Received 11 April 2017; Accepted 5 June 2017 Available online 06 June 2017 0921-4526/ © 2017 Elsevier B.V. All rights reserved.

low-temperature one demonstrates ferroelectric properties. Even at room temperature, there are two structurally inequivalent types of hydrogen bonds of different length, O2 − H2···O2 (R = 2.48 Å) and O3 − H3···O3 (R = 2.51 Å), phosphite HPO3 groups linked into chains along c-axis [13]. In both cases, the protons located there are disordered above TC over the two possible sites on hydrogen bonds whereas they become ordered in one of those sites below TC. This proton ordering in low-symmetry phase is a common feature of many hydrogen bonded ferroelectrics. However, there is an unique feature distinguishing hydrogen bonds of GPI among other ferroelectics of order-disorder type. As seen in Fig. 1, every O2 atom is involved only into a single O2 − H2···O2 hydrogen bond whereas every O3 atom is engaged into two different hydrogen bonds, i.e., O3 − H3···O3 and O3 − H6···N . This additional hydrogen bond may puzzle the mechanism of phase transition. Based on the lattice dynamics simulation of GPI in both structural phases [15] it was established that mean square displacements of O2 atoms in low-frequency range of 0–500 cm −1 are of one order of magnitude larger than those of O3 atoms. However, the H2 and H3 protons have almost similar thermal vibration amplitudes over the frequency range, 0–1000 cm −1. Hence, an attempt to shed light onto the GPI phase transition mechanism is the main motivation for the present study. In our former works, we succeeded in the efforts to properly describe the set of thermodynamic, elastic and piezoelectric properties of hydrogen-

Physica B 520 (2017) 164–173

I.R. Zachek et al.

Pseudo-spin variables

σq1

σq 4

describe the changes related to 2 → σqf reorientation of dipole moments of base units: dqf = μf 2 . Mean values σ

2

, …,

1

〈 2 〉 = 2 (na − nb ) imply the difference in occupancy of two possible proton positions on a hydrogen bond, na and nb. Hamiltonian of GPI proton may be presented as the sum of “seed” and pseudo-spin constituents. “Seed” energy Upseed effectively accounts for the crystal lattice and does not directly depend on proton configlshort urations. Pseudo-spin part of Hamiltonian considers short-range H lMF proton interactions near tetrahedral HPO3 groups and long-range H and also takes into account the interaction with cartesian components of electric field E1, E2 and E3. Therefore, Hamiltonian may be written as follows: l = NUseed + H lshort + H lMF − H

(2.1)

⎛ σq1 σq2 σq3 σq 4 ⎞ −E1 ∑ ⎜μ13x − μ24x − μ13x + μ24x ⎟− ⎝ 2 2 2 2 ⎠ q

Fig. 1. Crystal structure of GPI in paraelectric phase.

⎛ σq1 σq2 σq3 σq 4 ⎞ −E2 ∑ ⎜μ13y − μ24y + μ13y − μ24y ⎟− ⎝ 2 2 2 2 ⎠ q

bonded compounds with phase transitions of order-disorder type in the framework of pseudo-spin model for a proton subsystem within a twoparticle cluster approximation [16,17]. In the present work, we are aiming at elaborating a pseudo-spin model of GPI crystal capable of describing the set of thermodynamic properties observed experimentally. Constructing such a model we took into consideration the ground structural feature occurring at the phase transition, namely proton ordering within O2 − H2···O2 and O3 − H3···O3 hydrogen bonds. All other structural changes are treated by us effectively using adjustable model parameters. If our description adequately correlates with experiment, especially within the phase transition region, it should count in favor of order-disorder mechanism of the ferroelectric phase transition in GPI.

⎛ σq1 σq2 σq3 σq 4 ⎞ −E3 ∑ ⎜μ13z + μ24z − μ13z − μ24z ⎟. ⎝ 2 2 2 2 ⎠ q where v is unit cell volume, N is the total number of unit cells. In (2.1) μ13x, y, z , μ24x, y, z are effective dipole moments per one pseudo-spin; σqf is ycomponent of pseudo-spin operator, which corresponds to the proton located in the q-th cell on the f-th bond ( f = 1, 2, 3, 4 ). “Seed” Useed energy includes the elastic, piezolectric and dielectric constituents related to electric fields Ei(i = 1, 2, 3) and strains εi and εj (j = i + 3). cjjE 0(T ), eij0 , χiiε0 correspond to the so-called “seed” elastic constants, piezoelectric stresses and dielectric susceptibilities, respectively:

⎛ 1 Useed = v⎜⎜ ⎝2

2. Model Hamiltonian We consider the system of protons, based on O2 − H2···O2 and O3 − H3···O3 bonds, which form zigzag chains along crystallographic c→ axis. Dipole moments dqf (f = 1, 3) are ascribed to two H2 protons located on the neighboring chains along a-axis, whereas dipole → moments dqf (f = 2, 4) are related to two H3 protons on the same neighboring chains (see Fig. 2). In ferroelectric phase, the dipole → → → → moments compensate each other (dq1 with dq3, dq2 with dq4 ) along Z and X directions (X ⊥(b,c), Y ∥b, Z ∥c). However, along Y axis, they have non-zero components contributing to the total spontaneous polarization.

3

∑ i, i ′=1

1 ciiE 0(T )εiεi + ′ ′ 2

6

∑

3

cjjE 0(T )εj2 +

j =4

∑

E0 ciE5 0(T )εiε5 + c46 (T )

i =1

ε4ε6− 3 0 0 0 − ∑ e20iεiE2 − e25 ε5E2 − e140ε4E1 − e160ε6E1 − e34 ε4E3 − e36 ε6E3− i =1

(2.2)

⎞ 1 1 ε0 2 1 ε0 2 ε0 − χ11ε0 E12 − χ22 E2 − χ33 E3 − χ31 E3E1⎟. ⎠ 2 2 2 Hamiltonian of short-range interactions has the following form

Fig. 2. Orientations of vectors dqf in unit cell in ferroelectric phase.

165

Physica B 520 (2017) 164–173

I.R. Zachek et al.

⎞ ⎛σ σ σ σ ⎞⎛ lshort = − 2w ∑ ⎜ q1 q2 + q3 q 4 ⎟⎜δ R R + δ R + R , R ⎟ , H q q q c q ⎝ 2 2 2 2 ⎠⎝ ′ ′⎠ qq ′

⎛σ σ σ σ ⎞ lq(2)= −2w⎜ q1 q2 + q3 q 4 ⎟− H ⎝ 2 2 2 2 ⎠ y σq1 y1 σq1 y2 σq2 y σq 4 − − − 3 − 4 , β 2 β 2 β 2 β 2

(2.3)

where the first Kronecker symbol corresponds to interactions between the H2 or H3 protons near the tetrahedral HPO3 groups of type I in caxis running chains, the second – near the HPO3 groups of type II, Rc is the radius-vector of relative position of the hydrogen bond in the unit cell. Here, for the sake of simplicity, we introduced new notation of two types of HPO3 complexes. We refer to HPO3 tetrahedron along with two attached hydrogen bonds as I type tetrahedron, whereas the second type of HPO3 groups (II type) are accounted for here without O2 − H2···O2 and O3 − H3···O3 neighboring hydrogen bonds. Energy contributions of interactions between protons near the tetrahedra of different types, and mean values of the pseudo-spins <σqf >, which correspond to the tetrahedra of different type, are equal. Parameter w, which describes short-range interactions within -axis chains, is expanded linearly into a series with respect to strains εi, εj: 3 0

lqf(1) = − H

qq ′ ff ′

+

〈σqf 〉 〈σq f 〉 ′′ − 2 2

∑ Jff ′(qq′) qq ′ ff ′

1

1

1

1

1

1

1

1

β{Δ3 + 2 J11η3 + 2 J12η4 + 2 J13η1 + 2 J14η2−μ13x E1 + μ13y E2−μ13z E3},

y4=

β{Δ4 + 2 J22η4 + 2 J12η3 + 2 J24η2 + 2 J14η1 + μ24x E1−μ24y E2−μ24z E3},

(2.5)

∂Jff ∂εi

3

′ εi =

J ff0

l (2)

Sp σqf e−βHq

′

∑ ψff iεi + ∑ ψff jεj .

+

i =1

l (2)

Sp e−βHq

′

j =4

′

⎛

∑ ⎜⎝/ 1

σq1 2

(2.6)

η1 = 3

(2.7)

+ /2

σq2 2

+ /3

σq3 2

+ /4

η2 =

(2.9)

4

In Eq. (2.9), we used the following notations:

(2.10)

G = NUseed + NH 0− l (2)

3

6

∑ ln Spe−βHqf } − Nv ∑ σεi i − Nυ ∑ σjεj , f =1

. (2.16)

1 1 1 [sh (y1 + y2 + y3 + y4 ) ± sh (y1 + y2 − y3 − y4 ) D 2 2 1 1 2 2 +a sh (y1 − y2 + y3 − y4 ) ± a sh (y1 − y2 − y3 + y4 ) 2 2 1 1 ±ash (y1 − y2 + y3 + y4 ) ± ash (y1 + y2 + y3 − y4 ) 2 2

1 1 1 [sh (y1 + y2 + y3 + y4 ) ± sh (y1 + y2 − y3 − y4 )− D 2 2

(2.18)

1 1 +a 2ch (y1 − y2 + y3 − y4 ) + a 2ch (y1 − y2 − y3 + y4 )+ 2 2 1 1 +ach (y1 − y2 + y3 + y4 ) + ach (y1 + y2 + y3 − y4 )+ 2 2 1 1 +ach ( − y1 + y2 + y3 + y4 ) + ach (y1 + y2 − y3 + y4 ), 2 2

(2.11) l (1)

l (1)

Sp e−βHqf

y2 1 1 +ash ( − y1 + y2 + y3 + y4 ) + ash (y1 + y2 − y3 + y4 )] = th 4 , 2 2 2 1 1 D = ch (y1 + y2 + y3 + y4 ) + ch (y1 + y2 − y3 − y4 )+ (2.19) 2 2

To calculate the set of thermodynamic and dynamic characteristics of GPI, we utilized the two-particle cluster approximation. Within this approximation, the thermodynamic potential is given by:

−kBT ∑ {2ln Spe−βHq −

Sp σqf e−βHqf

1 1 −a 2sh (y1 − y2 + y3 − y4 ) ∓ a 2sh (y1 − y2 − y3 + y4 )∓ 2 2 1 1 ∓ash (y1 − y2 + y3 + y4 ) ± ash (y1 + y2 + y3 − y4 )+ 2 2

1 1 1 1 / 1= J11η1 + J12η2 + J13η3 + J14η4 , 2 2 2 2 1 1 1 1 / 2= J22η2 + J12η1 + J24η4 + J14η3, 2 2 2 2 1 1 1 1 / 3= J11η3 + J12η4 + J13η1 + J14η2 , 2 2 2 2 1 1 1 1 / 4= J22η4 + J12η3 + J24η2 + J14η1, 2 2 2 2

4

l (1)

=

y1 1 1 ±ash ( − y1 + y2 + y3 + y4 ) ± ash (y1 + y2 − y3 + y4 )] = th 3 , 2 2 2 (2.17)

(2.8)

σq 4 ⎞ ⎟. 2 ⎠

(2.15)

On the basis of Eq. (2.16), taking into account Eqs. (2.12) and (2.13), we obtain

6

1 1 1 1 J11(η12 + η32 ) + J22(η22 + η42 ) + J13η1η3 + J24η2η4+ 8 8 4 4 1 1 + J12(η1η2 + η3η4 ) + J14(η1η4 + η2η3), 4 4

q

1

y3=

where

q

1

Here, Δf are effective fields, created by the neighboring bonds from outside the cluster. In the cluster approximation, the fields Δf can be determined from the self-consistency condition, which states that the mean values of pseudo-spins 〈σqf 〉 calculated using two-particle and one-particle Gibbs distribution should coincide, that is,

〈σq f 〉 σqf ′′ . 2 2

lMF = NH 0 + H ls, H

ls = − H

1

β{Δ2 + 2 J22η2 + 2 J12η1 + 2 J24η4 + 2 J14η3−μ24x E1−μ24y E2 + μ24z E3},

yf = βΔf + yf .

Finally, one may write the Hamiltonian (2.5) as

H 0=

1

y2=

(2.14)

(2.4)

j =4

∑ Jff ′(qq′)

′

1

1

1

1

Fourier transforms of interaction constants at q = 0 J ff ′ = ∑q ′ J ff ′(qq′) may be linearly expanded into a series with respect to strains εi, εj:

Jff = ′

(2.13)

y1= β{Δ1 + 2 J11η1 + 2 J12η2 + 2 J13η3 + 2 J14η4 + μ13x E1 + μ13y E2 + μ13z E3},

lMF is the mean field Hamiltonian of the long-range dipol-dipol H interactions and indirect (through lattice vibrations) proton interactions:

J ff0

,

6

i =1

lMF = 1 H 2

β 2

and the following notation was used in (2.12) and (2.13):

∑ δiεi + ∑ δjεi .

w=w +

yf σqf

(2.12)

i =1

⎛ − 1 ⎜w0 +∑3

j =4

δiεi +∑6

⎞ δjεj ⎟

i =1 j =4 ⎠. where a = e kBT ⎝ Combining Eqs. (2.14), (2.15), (2.17) and (2.18), one may exclude parameters Δf and come finally to relations:

lq(2), H lqf(1) are two-particle and one-particle Hamiltonians, rewhere H spectively: 166

Physica B 520 (2017) 164–173

I.R. Zachek et al.

y1=

1 + η1 β 1 ln + βν11η1 + βν12η2 + βν13η3 + βν14η4 + (μ13x E1 + μ13y E2 2 1 − η1 2

From equilibrium conditions

1 ⎛ ∂g ⎞ 1 ⎛ ∂g ⎞ ⎜ ⎟ = 0, ⎜⎜ ⎟⎟ v ⎝ ∂εi ⎠ v ⎝ ∂εj ⎠

+ μ13z E3), y2= βν12η1 +

Ei

1 + η2

1 ln + βν22η2 + βν14η3 + βν24η4 2 1 − η2

we obtain equations for strains εi, εj

β ( − μ24x E1 − μ24y E2 + μ24z E3), 2 1 + η3 β 1 y3= βν13η1 + βν14η2 + ln + βν11η3 + βν12η4 + ( − μ13x E1 + μ13y E2 2 1 − η3 2 +

2δl 2δ + l Mε υ vD ψ ψ ψ ψ − 11 l (η12 + η32 )− 13 l η1η3− 22 l (η22 + η42 )− 24 l η2η4 8v 4v 8v 4v ψ ψ − 12 l (η1η2 + η3η4 )− 14 l (η1η4 + η2η3), (l = 1, 2, 3, 5) 4v 4v 2δ 2δ E0 E0 0 σ4= c44 ε4 + c46 ε6−e140E1−e34 E3− 4 + 4 Mε υ vD ψ ψ ψ ψ − 114 (η12 + η32 )− 134 η1η3− 224 (η22 + η42 )− 244 η2η4 8v 4v 8v 4v ψ ψ − 124 (η1η2 + η3η4 )− 144 (η1η4 + η2η3), 4v 4v 2δ 2δ E0 E0 0 0 σ6= c46 ε4 + c66 ε6−e16E1−e36E3− 6 + 6 Mε υ vD ψ136 ψ226 2 ψ ψ116 2 2 2 ηη− − (η + η3 )− (η + η4 )− 246 η2η4 8v 1 4v 1 3 8v 2 4v ψ126 ψ146 − (η η + η3η4 )− (η η + η2η3). 4v 1 2 4v 1 4 σl = clE1 0ε1 + clE2 0ε2 + clE3 0ε3 + clE5 0ε5−e20l E2−

− μ13z E3), y4= βν14η1 + βν24η2 + βν12η3 + +

1 + η4 1 ln + βν22η4 2 1 − η4

β x (μ E1 − μ24y E2 − μ24z E3), 2 24 J ff

where νff = 4 ′ . ′ In paraelectric phase, at the absence of external electric fields and mechanical stresses, the mean values of pseudo-spins η1 = η2 = η3 = η4 = 0 . In ferroelectric phase, at zero electric fields Ei = 0 and mechanical stresses: σj = 0, η1 = η3 = η13, η2 = η4 = η24 ,

1 [sh(y13 + y24 ) + a 2sh(y13 − y24 ) + 2ashy13], D͠ 1 η24= [sh(y13 + y24 ) − a 2sh(y13 − y24 ) + 2ashy24] D͠ η13=

1 1 Mε= 2a 2ch (y1 − y2 + y3 − y4 ) + 2a 2ch (y1 − y2 − y3 + y4 )+ 2 2 1 1 +ach (y1 − y2 + y3 + y4 ) + ach (y1 + y2 + y3 − y4 )+ 2 2 1 1 +ach ( − y1 + y2 + y3 + y4 ) + ach (y1 + y2 − y3 + y4 ). 2 2

where the following notations are used:

1 + η13 1 + η24 1 1 ln + βν1+η13 + βν2+η24 , y24 = βν2+η13 + ln 2 1 − η13 2 1 − η24 + βν3+η24 ,

⎛ 3 ν1+= ν10+ + ⎜⎜∑ ψ1+i εi + ⎝ i =1

ν2+=

ν20+

⎛ 3 + ⎜⎜∑ ψ2+i εi + ⎝ i =1

∑ j =4

6

∑ j =4

= ⎛ ν3+= ν30+ + ⎜⎜∑ ψ3+i εi + ⎝ i =1 3

6

∑ j =4

1 ⎛ ∂g ⎞ Using the definition of polarization, Pi = − v ⎜ ∂E ⎟, we come to the ⎝ i⎠ components of polarization

⎞ 1 ψ1+j εj ⎟⎟ , ν10+ = (J110 + J130); ψ1+i 4 ⎠ 1 = (ψ11i + ψ13i ), 4 ⎞ 1 ψ2+j εj ⎟⎟ , ν20+ = (J120 + J140 ); ψ2+i 4 ⎠

6

1

+ 2v [μ13y (η1 + η3)−μ24y (η2 + η4 )], 0 0 ε0 ε0 P3 = e34 ε4 + e66 ε6 + χ33 E3 + χ31 E1 +

j

i =1

+

β − {(μ x )2[D͠ λ 24 − (λ13λ 24 − λ 2 )φ24 ]+ 2υΔ1,3 13

−2μ13x μ24x [D͠ λ + (λ13λ 24 − λ 2 )βν2−]},

δjεi )+

j =4

2 − Δ1,3 = D͠ − D͠ (λ 24φ13− + λ13φ24 + 2λβν2−)+

4 3 ⎞ ⎛ 1 +2kBT ln2 − kBT ∑ ln ⎜1 − ηf2⎟ − 2kBT ln D − Nv ∑ σε i i − Nυ ⎠ ⎝ 2 f =1 i =1

− +(λ13λ 24 − λ 2 )[φ13−φ24 − (βν2−)2 ],

6

∑

(3.4)

+(μ24x )2[D͠ λ13 − (λ13λ 24 − λ 2 )φ13−]−

6

∑

1 z [μ (η − η3) + μ24z (η2 − η4 )]. 2v 13 1

⎛ ∂P ⎞ χ11Tε (0) = lim ⎜ 1 ⎟ = χ11ε0 + E1→0⎝ ∂E1 ⎠ ε

To calculate the dielectric, piezoelectric and elastic characteristics of GPI, we use thermodynamic potential per unit cell obtained in the twoparticle cluster approximation:

δiεi +

(3.3)

For simplicity we restrict our analysis to the case of zero mechanical stress and zero electric field. Components of static isothermal dielectric susceptibility of a mechanically clamped GPI crystal are given by:

3. Thermodynamic characteristics of GPI

∑

0 0 e21 ε1 + e22 ε2 +

P2 =

1 (ψ + ψ14i ), 4 12i

3

1 [μ x (η −η3)−μ24x (η2−η4 )], 2v 13 1 ε0 0 0 e23 ε3 + e25 ε5 + χ22 E2

ε0 P1 = e140ε4 + e160ε6 + χ11ε0 E1 + χ31 E3 +

⎞ 1 0 0 ψ3+j εj ⎟⎟ , ν30+ = (J22 + J24 ); ψ3+i 4 ⎠ 1 = (ψ22i + ψ24i ). 4

G = Useed + H 0 − 2(w0 + g= N

(3.2)

where

(2.20)

D͠ = ch(y13 + y24 ) + a 2ch(y13 − y24 ) + 2achy13 + 2achy24 + a 2 + 1,

y13=

=0 Ei, σi

⎛ ∂P ⎞ Tε ε0 (0) = lim ⎜ 2 ⎟ = χ22 χ22 + E2→0⎝ ∂E2 ⎠ ε

σjεj .

j =4

j

(3.1) 167

(3.5)

Physica B 520 (2017) 164–173

I.R. Zachek et al.

β + {(μ13y )2 [D͠ ϰ13 − (ϰ13ϰ24 − ϰ2)φ24 ]+ 2υΔ2 +(μ y )2 [D͠ ϰ24 − (ϰ13ϰ24 − ϰ2)φ +]−

β [(ψ η + ψ2iη24 )τ1ψ + (ψ2iη13 + ψ3iη24 )τ2ψ − 2δiτ1δ ], Δ2 1i 13 β γ24i= [(ψ η + ψ2iη24 )τ2ψ + (ψ2iη13 + ψ3iη24 )τ3ψ − 2δiτ2δ ], Δ2 1i 13 + τ1ψ = D͠ ϰ13 − (ϰ13ϰ24 − ϰ2)φ24 , τ2ψ = D͠ ϰ + (ϰ13ϰ24 − ϰ2)βν2+,

+

24

γ13i=

13

−2μ13y μ24y [D͠ ϰ + (ϰ13ϰ24 − ϰ2)βν2+]} 2 + Δ2 = D͠ − D͠ (ϰ13φ13+ + ϰ24φ24 + 2ϰβν2+)+

τ3ψ = D͠ ϰ24 − (ϰ13ϰ24 − ϰ2)φ13+,

+ +(ϰ13ϰ24 − ϰ2)[φ13+φ24 − (βν2+)2 ],

Tε χ33 (0)

+ + τ1δ = [D͠ − (ϰ24φ24 + ϰβν2+)]ρ13 + (ϰφ24 + ϰ13βν2+)ρ24 ,

τ2δ = [D͠ − (ϰ13φ13+ + ϰβν2+]ρ24 + (ϰφ13+ + ϰ24βν2+)ρ24 ,

⎛ ∂P ⎞ ε0 = lim ⎜ 3 ⎟ = χ33 + E3→0⎝ ∂E3 ⎠ ε

ρ13= [a 2sh(y13 − y24 ) + ashy13] − η13M , (3.6)

j

ρ24 = [ − a 2sh(y13 − y24 ) + ashy24] − η24M , M = a 2ch(y13 − y24 ) + achy13 + achy24 + a 2.

β − + {(μ z )2 [D͠ λ 24 − (λ13λ 24 − λ 2 )φ24 ]+ 2υΔ1,3 13 +(μ24z )2 [D͠ λ13 +2μ13z μ24z [D͠ λ

− (λ13λ 24 − + (λ13λ 24 −

λ )φ13−]+ λ 2 )βν2−]},

⎛ ∂P ⎞ ε ε0 χ31 = lim ⎜ 3 ⎟ = χ31 + E1→0⎝ ∂E1 ⎠ ε

⎛ ∂E ⎞ e e h 2is = − ⎜ 2 ⎟ = 2εis , h 25 s = 25ε s . χ22 s χ22 s ⎝ ∂εi ⎠P 2

⎛ ∂σ ⎞ cijE = ⎜ i ⎟ = cijE0− ⎝ ∂εi ⎠E

β − {μ x μ z [D͠ λ 24 − (λ13λ 24 − λ 2 )φ24 ]− 2υΔ1,3 13 13

2β {(ψ1iη13 + ψ2iη24 )(ψ1jη13 + ψ2jη24 )τ1ψ + vΔ2 +[(ψ1iη13 + ψ2iη24 )(ψ2jη13 + ψ3jη24 ) + (ψ2iη13 + ψ3iη24 )(ψ1jη13 + ψ2jη24 )]τ2ψ +

−(μ24x μ13z − μ13x μ24z )[D͠ λ + (λ13λ 24 − λ 2 )βν2−]}.

−

with the following notations:

+(ψ2iη13 + ψ3iη24 )(ψ2jη13 + ψ3jη24 )τ3ψ }+

1 1 ± + βν1±, φ24 = + βν3±, 2 1 − η132 1 − η24

⎛ 3 ν1−= ν10− + ⎜⎜∑ ψ1−i εi − ⎝ i =1

∑

⎛ 3 ν2−= ν20− + ⎜⎜∑ ψ2−i εi − ⎝ i =1

∑

6 j =4

6 j =4

4βδi [(ψ η + ψ2jη24 )τ1δ + (ψ2jη13 + ψ3jη24 )τ2δ ]+ υΔ2 1j 13 4βδj + [(ψ η + ψ2iη24 )τ1δ + (ψ2iη13 + ψ3iη24 )τ2δ ]− υΔ2 1i 13 8βδiδj + − {(ρ φ + + ρ24 βν2+)τ1δ + (ρ24 φ24 + ρ13βν2+)τ2δ}− υD͠ Δ2 13 13 +

⎞ 1 ψ1−j εj ⎟⎟ , ν10− = (J110 − J130); ψ1−i 4 ⎠ 1 = (ψ11i − ψ13i ), 4 ⎞ 1 ψ2−j εj ⎟⎟ , ν20− = (J120 − J140 ); ψ2−i 4 ⎠ =

−

6

∑

2 ͠ D, ϰ24= ch(y13 + y24 ) + a 2ch(y13 − y24 ) + 2achy24 − η24

d 2i =

2

c12E c13E c15E ⎞ ⎟ E E E⎟ c22 c23 c25 ⎟, E E E c23 c33 c35 ⎟ E E E⎟ c25 c35 c55 ⎠

m S E = (m C E )−1,

∑ sijEe2j , (i, j = 1, 2, 3, 5), j

Using Eqs. (3.2) and (3.3) we obtain the expressions for isothermal coefficients of piezoelectric stress e2j:

⎛ ∂P ⎞ μy μy 0 + 13 γ135 − 24 γ245, e25 = ⎜ 2 ⎟ = e25 v v ⎝ ∂ε5 ⎠E

{[2a 2ch(y13 − y24 ) + achy13 + achy24 + 2a 2]D͠ − 2 M2. }

isothermal constants of piezoelectric strain

ϰ= ch(y13 + y24 ) − a 2ch(y13 − y24 ) − η13η24D͠ .

2

2 υD͠

⎛c E ⎜ 11 ⎜c E m E C = ⎜ 12E ⎜ c13 ⎜ E ⎝ c15

ϰ13= ch(y13 + y24 ) + a 2ch(y13 − y24 ) + 2achy13 − η132D͠ ,

⎛ ∂P ⎞ μy μy e2i = ⎜ 2 ⎟ = e20i + 13 γ13i − 24 γ24i, v v ⎝ ∂εi ⎠E

4βδiδj

Other dielectric, piezoelectric and elastic characteristics of GPI can be found using the expressions established above. In particular, the matrix of isothermal elastic compliance at a constant field sijE, which is reciprocal to the matrix of elastic constants cijE is as follows:

1 (ψ − ψ14i ), 4 12i

⎞ 1 0 0 ψ3−j εj ⎟⎟ , ν30− = (J22 − J24 ); ψ3−i 4 j =4 ⎠ 1 = (ψ22i − ψ24i ). 4 λ13= 1 + a 2 + 2achy13, λ 24 = 1 + a 2 + 2achy24 , λ = 1 − a 2 , ⎛ 3 ν3−= ν30− + ⎜⎜∑ ψ3−i εi − ⎝ i =1

(3.11)

2

−μ24x μ24z [D͠ λ13 − (λ13λ 24 − λ 2 )φ13−]−

φ13±=

(3.10)

Proton contribution to elastic constants of GPI is found by differentiating (3.2) over strains at a constant field:

(3.7)

j

+

Constants of piezoelectric stress are obtained by differentiation of the electric field, found from Eq. (3.3), over strains at constant polarization:

2

(3.12)

isothermal dielectric susceptibility of mechanically free crystal σ ε χ22 = χ22 +

∑ e2id2i + e25d25, i

(3.8)

(3.13)

the constants of piezoelectric strain

g2i = (3.9)

d 25 d 2i σ , g25 = σ . χ22 χ22

(3.14)

Molar entropy of the proton subsystem (here, R is the gas constant):

where 168

Physica B 520 (2017) 164–173

I.R. Zachek et al.

Tc, K

350

“seed” coefficients of piezoelectric stress eij0 ; “seed” elastic constants cijE0 . To establish these parameters, we used the experimental temperature dependence of the following physical characteristics: Ps(T) [18], σ Cp(T) [6], ε11σ , ε33 [19], d21, d23 [20]. The hydrostatic pressure dependence of phase transition temperature Tc(p) [21] was used as well. Static dielectric permittivity was analyzed by us for a set of partially deuterated crystals GPI 1− x DGPIx with different deuterium content. Fully deuterated analogue of GPI has the following formulae ND3CH2COOD·D2PO3. It means that the hydrogen atoms are substituted by deuterium in O2 − D2 ···O2 and O3 − D3···O3 hydrogen bonds which play a key role in the ferroelectric phase transition. Deuteron content x dependence of phase transition temperature Tc of GPI 1− x DGPIx was published in two papers [22,23]. The agreement between these data is rather poor, especially for high deuterium concentration, x≥ 0.3 (see Fig. 3). Since the papers [18,22] present the concentration dependence of spontaneous polarization and dielectric permittivity of the mixed GPI 1− x DGPIx compound we decided to utilize in our further analysis the content x dependence published in [22]. As seen in Fig. 3, these experimental data [22] may be nicely approximated by the following relation Tc(x ) = 225(1 + 0.382x + 0.193x 2 )K (blue line). Hereinafter, this relation is accepted by us as the concentration law for Tc in GPI 1− x DGPIx. We admitted the linear dependence on deuterium x concentration of the unit cell volume of a mixed compound, i.e., υ(x ) = υH (1 − x ) + υDx , where υH=0.6015 ·10−21 cm3 [13], υD=0.6139 ·10−21 cm3 [23] are the unit cell volumes of fully protonated and deuterated compound, respectively. The analysis of numerical data shows that thermodynamic characteristics mostly depend on two linear combinations of long-ranged interactions, i.e., ν 0+ = ν10+ + 2ν20+ + ν30+ and ν 0− = ν10− + 2ν20− + ν30− being practically independent (deviation < 0.1%) of the particular value of each νf0± at the given ν0+ν0-. The same condition regarding the relevance of linear combinations, ψi+ = ψ1+i + 2ψ2+i + ψ3+i and ψi− = ψ1−i + 2ψ2−i + ψ3−i , is valid for ψ fi± parameters. We have got the

300

250

200

0

0.2

0.4

0.6

0.8

x

Fig. 3. Phase transition temperature Tc vs deuterium concentration x obtained in GPI DGPIx: [22]; [23].

1− x

S=−

R ⎛ ∂g ⎞ R ⎜ ⎟ = {−2ln2 + ln(1 − η13) + ln(1 − η24 ) + 2lnD͠ − 4 ⎝ ∂T ⎠η, ε 4 i

(3.15)

−2(βν1+η13 + βν2+η24 )η13 − 2(βν2+η13 + βν3+η24 )η24 +

4w ⎫ M ⎬, TD͠ ⎭

The molar heat capacity of the proton subsystem of GPI crystals can be found from the entropy Eq. (3.15):

⎛ ∂S ⎞ ΔC σ = T ⎜ ⎟ . ⎝ ∂T ⎠σ

(3.16)

4. Comparison between theory and experiment. Discussion

optimal values of these combinations ν 0+ / kB = 10.57 K, ν 0− / kB =−0.8 K with the following values of νf0± , i.e., ν∼10+ = ν∼20+ = ν∼30+ = 2.643 K, ∼+ = ∼+ = 87.9 K, ψ ν∼ 0− = ν∼ 0− = ν∼ 0−=0.2 K, where ν∼0± = ν 0± / k and ψ

To calculate the temperature evolution of dielectric, elastic and piezoelectric characteristics of GPI, we need to set certain values of the following parameters: parameter of short-range interactions w0; parameters of long-range interactions νf0± (f = 1, 2, 3); deformation potentials δi, δ5, ψ11i, ψ13i, ψ12i, ψ14i, ψ22i, ψ24i; effective dipole moments μ13x ; μ24x ; μ13y ; μ24y ; μ13z ; μ24z ; “seed” dielectric susceptibilities χijε0 ;

Fig. 4. Temperature dependence of spontaneous polarization PSGPI (a) – 0.47,

[18], 5 - 0.67,

[18], 6 - 0.75,

[18],

[19],

1

2

f

3

f

B

f1

f2

∼+ = 149.1 K, ψ ∼+ = 21.3 K, ψ ∼+ = 143.8 K, ψ ∼+ = 103.8 K, ψ ∼− 237.0 K, ψ f4 f5 f6 f3 fi ∼± = ψ ± / k . The parameters ν 0± are independent of the = 0 K, where ψ fi

fi

B

f

deuterium x content, i.e., νf0±(x ) = const, whereas the ψ fi± parameters linearly decrease with x increasing according to the following law ψ fi±(x ) = ψ fi±(0)(1 − 0.82x ). It appears that η13 > η24 at ν10+ > ν30+ and

[20],

[18]. The lines correspond to the calculated data.

169

[1] and GPI 1− x DGPIx (b)– 1 - x=0.00,

[18], 2 - 0.16,

[18], 3 - 0.31,

[18], 4 -

Physica B 520 (2017) 164–173

I.R. Zachek et al.

Fig. 5. Temperature dependence of static dielectric permittivities of GPI (a): (b): 1 - x=0.00,

[22], 2 - 0.16,

[22], 3 - 0.31,

[22], 4 - 0.47,

- 1 kHz [19],

[22], 5 - 0.67,

- 10 kHz [22],

[22], 6 - 0.75,

- 10 kHz [21],

[22]. The lines correspond to the calculated data.

η13 < η24 at ν10+ < ν30+ . Since ferroelectric phase transition is a transition of the second order, based on the extremum condition of the reciprocal dielectric susceptibility (3.5), one may readily obtain the equation for phase transition temperature, Δ2p (Tc ) = 0 , which allows us to find the w0 parameter for various deuterium concentrations x. Finally, the w0(x ) dependence may be approximated with high accuracy (< 0.1%) by means of cuspidal function, w0(x )/ kB = 820(1 + 0.529x + 0.286x 2 )K . The deformation potentials, δj, are the coefficients of w0 expansion into a series with respect to strain εj (see Eq. (2.4)), and have the ∼ ∼ ∼ ∼ ∼ following values, δ 1 = 500 K, δ 2 = 600 K, δ 3 = 500 K, δ 4 = 150 K, δ 5 = ∼ ∼ 100 K, δ 6 = 150 K; δ i = δi / kB . There is a linear concentration x ∼ ∼ dependence of all δi parameters, δi(x ) = δi(0)(1 − 0.82x ). In paraelectric phase, the effective dipole moments of 1 (3) and 2 (4) protons are independent of the deuterium content and adopt certain values along the Cartesian axes, → μ 13 = (0.4, 4.02, 4.3)·10−18 −18 → esu cm, μ 24 = (−2.3, −3.0, 2.2)·10 esu cm, respectively. It appears, in ferroelectric phase, that μ13y component grows with an increasing x according to the linear relation, i.e. μ13y ferro (x ) = 3.82(1 + 0.062x )·10−18 esu cm. Using experimental structural data in both temperature phases [12], we calculated the displacements of protons labeled as 1 and 2 (see

Fig. 6. Temperature dependence of transverse permittivities ε11, ε33 and ε31 of GPI; [19], [19]. The lines correspond to the calculated data.

ε 0 −1 σ 0 −1 − ε22 ) - 1 and (ε22 − ε22 ) - 2 of GPI (a): Fig. 7. Temperature dependence of reciprocal dielectric permittivity (ε22

(b)– 1 - x=0.00,

[22], 2 - 0.31,

[22], 3 - 0.47,

[22], 4 - 0.67,

[22], 5 - 0.75,

ε σ (T ) and GPI 1− x DGPIx - 1 MHz [20], 1 - ε22 (T ) , 2 - ε22

[22], 6 - 0.75,

170

[22],

[20],

[19],

ε 0 −1 − ε22 ) of GPI 1− x DGPIx [21] and (ε22

[22]. The lines correspond to the calculated data.

Physica B 520 (2017) 164–173

I.R. Zachek et al.

Now let us discuss the obtained results. The comparison between the calculated temperature dependence of spontaneous polarization PS of protonated GPI and experimental data [18–20,1] is presented in Fig. 4a and a similar dependence for the mixed GPI 1− x DGPIx compound with different deuteron concentration x is shown in Fig. 4b. As seen in Fig. 4a, there is a disagreement between the experimental data published in different papers. The saturation polarization increases linearly with an increasing deuteron concentration x. Fig. 5a shows the calculated temperature dependence of static ε (T ) and free dielectric permittivity of the mechanically clamped ε22 σ ε22(T ) GPI crystal along with the experimental data [22,1,20,21,19], while Fig. 5b demonstrates the corresponding theoretical and experimental data [22] for GPI 1− x DGPIx crystals. As one can see from Fig. 5, ε ε (T ) = 1 + 4πχ22 ) well in paraelectric phase, the results of calculation (ε22 quantitatively agree with experimental data [22,20,21,19] for different deuteron concentrations x. However, in ferroelectric phase, the experimental data have much higher values compared to the calculated ones. This is due to a significant domain-wall motion contribution to the dielectric properties at low-frequencies of external electric field. As seen in Fig. 5, the domain-wall contribution to ε22 does increase with a decreasing frequency of the electric field. However, the low-frequency domain-wall dynamics is out of our current research. The temperature dependence of transverse dielectric permittivities ε11, ε33 and ε31 is presented in Fig. 6. Our theoretical results well describe the temperature dependence of ε11 permittivity in the whole temperature range investigated and of ε33 permittivity only in ferroelectric phase. There is a significant deviation of the calculated ε33(T) curve from the corresponding experimental data in paraelectric phase. The reason of this issue may be found in the search of physical origin of this unexpectedly high room temperature experimental value of ε33, i.e., ∼190. Normally, the permittivity of a dielectric material in non-polar phase does not exceed 10. In our opinion, such a large value observed along the c (Z) direction is related to the GPI crystal structure. Since both the phosphite chains and the hydrogen bonds are directed approximately along the c axis (see Fig. 2), there is an enhanced correlation of long-ranged interionic interactions exactly along the c direction. It causes a significant space anisotropy of polarizability of both ionic complexes and hydrogen bonds which results in high effective polarization exactly along c axis. Since our treatment is basically focused on the proton dynamics and effectively accounts all other interatomic interactions it may explain a significant disagreement between the calculated and experimental values of ε33 at room temperature (Fig. 6).

−1 −1 , ε22 and Fig. 8. Temperature dependence of the reciprocal transverse permittivities ε11 −1 ε33 of GPI; [19], [19], [22]. The lines correspond to the calculated data.

Fig. 2) to the equilibrium ferroelectric coordinates relative to the r1 = (−0.016, −0.495, centres of hydrogen bonds, as follows Δ→ r2 = (0.389, 0.383, −0.147)Å. As seen, there is a −0.160)Å, Δ→ quantitative correlation between our model parameters, → μ 13, → μ 24, r1 , Δ→ r2 . The absence of a quantitaand experimental displacements, Δ→ tive correspondence between our model and experimental data may be due to the fact that model dipole moments effectively account not only for the proton displacements but also for the distortions of phosphite HPO3 and glycine NH3CH2COOH complexes appearing at a phase transition. We came to the following “seed” coefficients of piezoelectric stress, dielectric susceptibilities and elastic constants. esu ε0 ε0 ε0 ε0 0 0 0 0 = 0.5, χ31 = 0.0 ; e21 = e22 = e23 = e25 = 0.0 2 ; χ11 = 0.1, χ22 = 0.403, χ33 cm

c110E = 26.91·1010

dyn

cm 10 dyn

= 3.91·10

2,

c12E 0 = 14.5·1010

, c13E 0 = 11.64·1010

dyn cm 2

, c15E 0 ,

cm 2

E0 = (64.99 − 0.04(T − Tc ))·1010 c22

= 5.64·1010 E0 c33

dyn cm 2

E0 , c23 = 20.38·1010

dyn cm 2

E0 , c25 ,

dyn

cm 2 10 dyn

= 24.41·10

cm 2 10 dyn

E0 c44 = 15.31·10

dyn cm 2

cm 2

dyn dyn E0 , c55 = 8.54·1010 2 , cm 2 cm dyn dyn E0 1.1·1010 2 , c66 = 11.88·1010 2 . cm cm

E0 , c35 = − 2.84·1010 E0 , c46 =−

Fig. 9. Temperature dependence of the coefficients of piezoelectric strain d21–1, stress e21–1, e22–2, e23–3, e25–4 (b) of GPI. Lines are the theoretical results.

[20], d22–2, d23–3,

171

[20], d25–4 (a) and temperature dependence of the coefficients of piezoelectric

Physica B 520 (2017) 164–173

I.R. Zachek et al.

Fig. 10. Temperature dependence of the constants of piezoelectric stress h21–1, h22–2, h23–3, h25–4 (a) and of the constants of piezoelectric strain g21–1, g22–2, g23–3, g25–4 (b) of GPI.

deuteron x concentrations. In ferroelectric phase, the calculated σ (0, T ))−1 is half the permittivity in paraelectric phase, permittivity (ε22 which corresponds to the Curie-Weiss law. Fig. 8 demonstrates the temperature dependence of reciprocal −1 −1 −1 , ε22 and ε33 of GPI crystal in a wide dielectric permittivities ε11 temperature range.

Fig. 7a shows the temperature dependence of reciprocal dielectric ε 0 −1 (0, T ) − ε22 ) and mechanically permittivity of mechanically free (ε22 σ 0 −1 clamped (ε22 (0, T ) − ε22 ) protonated GPI, while Fig. 7b - of the mixed GPI 1− x DGPIx crystals. As seen from Fig. 7a, in paraelectric phase, the σ (0, T ))−1 of GPI generally well agree with results of calculations for (ε22 experimental data [22,1,20,21,19] and with data [22] at different

Fig. 11. Temperature dependence of elastic constants of GPI crystal: c11 − 1, c12 − 2, c13 − 3, c22 − 4, c23 − 5, c33 − 6, c15 − 7, c25 − 8, c35 − 9, c55 − 10, c44 ;

– experimental data [5];

− 11, c46 − 12, c66 − 13 – c22 = ρv22 , [4];

– c22 = ρv22 [3] (a); temperature dependence of elastic compliances sij of GPI: s11 − 1, s12 − 2, s13 − 3, s22 − 4, s23 − 5, s33 − 6, s15 − 7, s25 − 8, s35 − 9, s55 − 10, s44

− 11, s46 − 12, s66 − 13 (b). Lines correspond to calculated results.

172

Physica B 520 (2017) 164–173

I.R. Zachek et al.

Cp(T ) − ΔCp(T ) near the phase transition point, where ΔCp(T ) is the quasi-spin contribution to the heat capacity near TC. The lattice contribution Clattice in the vicinity of Tc is well approximated by the straight line Clattice = C0 + C1(T − Tc ), where C0 = 147.5J /(mol · K ), C1 = 0.46J /(mol · K 2 ). Our calculation well agrees with the experimental temperature dependence of Cp(T) [6]. 5. Conclusion In our research we presented the modified model of GPI crystal based on two-particle cluster approximation which properly describes a rich set of various experimental thermodynamic characteristics. Crystal structure of GPI is very complicated consisting of various kinds of crystal interactions. Our model is based on the most characteristic feature of the GPI, namely the proton ordering within O2 − H2···O2 and O3 − H3···O3 hydrogen bonds occurring at structural phase transition. All other interatomic interactions and structural changes accompanying this ferroelectric transition were effectively accounted for using adjustable parameters. However, such proton ordering is also accompanied by intricate internal deformations of both phosphite HPO3 and glycinium groups NH3CH2COOH. Consequently, generally speaking, the quasi-spin formalism used in our model should be applicable not only to the proton dynamics in hydrogen bonds but also seemingly to all structural units of GPI undergoing the order-disorder-like dynamics at structural phase transition. We succeeded in describing the temperature dependence of dielectric, elastic, piezoelectric and thermodynamic properties. For most of them, we have got very good agreement between the model treatment and experiment. For some of them, e.g., for ε33(T) in paraelectric phase, the agreement was worse. In our opinion, the disagreement between our theoretical results and experimental data clearly indicates the limit of the model presented here and shows hopefully the eventual ways for its improvement.

Fig. 12. Temperature dependence of heat capacity of GPI, – [6]. Solid line presents the result of calculation, dashed line corresponds to lattice contribution.

The temperature dependence of the coefficients of piezoelectric strains d2i and d25 and the coefficients of piezoelectric stresses e2i and e25 is presented in Fig. 9, while the temperature dependence of the constants of piezoelectric stresses h2i, h25 and of piezoelectric strains g2i, g25 – in Fig. 10. In paraelectric phase, these coefficients should be equal to zero according to the inversion symmetry of crystal above TC. In ferroelectric phase, the values of d2i, d25 and e2i, e25 increase with an increasing temperature and reach a maximum at Curie point Tc. Temperature dependence of piezoelectric moduli h2i, h25 and g2i, g25 manifests a very broad maximum in ferroelectric phase centered near ∼220 K. As seen in Fig. 9, experimental data for d21 and d23 [20] demonstrate a rather unexpected behavior with broad diffusive peak near Curie temperature spreading its nonzero wing far above the phase transition Tc temperature. It is worth noting that according to the same paper [20], spontaneous polarization was not equal to zero above Tc thus indicating either some peculiarities of the GPI sample or some specific features of the experimental setup used in paper [20]. Fig. 11a shows the temperature dependence of elastic moduli cij at a constant electric field. All elastic moduli are temperature independent in paraelectric phase and undergo certain jumps down at Curie point TC showing a smooth increase with a decreasing temperature in the ferroelectric phase. In the same Figure, we present the experimental room temperature data of [5]. Since the papers [4,3] reported a temperature evolution of ultrasonic wave υ2 velocity, we recalculated the elastic modulus c22 using the relation c22 = ρυ22 , where the crystal density ρ = 1, 722 g / cm3 is taken from [5]. As seen in this Figure, there is a good agreement between the theory and experiment at room temperature for all elastic moduli. There is also a quite acceptable agreement between our calculation and experiment [4] for thermal evolution of c22 modulus with the exception of the close vicinity of phase transition point where the jump down of c22 is detected in our research and is not visible in experiment [4]. A quite discernible discontinuity down of c22 near Tc was clearly measured in [3] showing also the temperature evolution of c22 compared to our data giving, however, the room temperature value of c22 significantly lower than that published in [4,5]. The temperature dependence of sij elastic compliances of GPI crystal is given in Fig. 11b. Temperature dependence of GPI heat capacity is shown in Fig. 12. Dashed line in this Figure demonstrates the effective lattice Clattice contribution to Cp estimated by us as the average difference

References [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15] [16] [17] [18] [19] [20] [21] [22] [23]

173

R. Tchukvinskyi, R. Cach, Z. Czapla, S. Dacko, Phys. Stat. Sol. (a) 165 (1998) 309. J. Baran, G. Bator, R. Jakubas, M. Sledz, J. Phys.: Condens Matter 8 (1996) 10647. J. Furtak, Z. Czapla, A.V. Kityk, Z. für Nat. A 52 (1997) 778. K. Lapsa, M. Drozdowski, P. Ziobrowski, L. Szczepanska, Ferroelectrics 239 (2000) 87. A. Deepthy, H.L. Bhat, A.V. Alex, J. Philip, Phys. Rev. B 62 (2000) 8752. F. Shikanai, J. Hatori, M. Komukae, Z. Czapla, T. Osaka, J. Phys. Soc. Jpn. 73 (2004) 1812. G. Bhoopathi, V. Jayaramakrishnan, K. Ravikumar, T. Prasanyaa, S. Karthikeyan, Mater. Sci. Pol. 31 (2013) 1. T. Runka, M. Kozelski, M. Drozdowski, L. Szczepanska, Ferroelectrics 239 (2000) 995. M. Sledz, J. Baran, J. Mol. Struct. 706 (2004) 15. I. Stasyuk, Z. Czapla, S. Dacko, O. Velychko, J. Phys.: Condens Matter 16 (2004) 1963. I. Stasyuk, O. Velychko, Ferroelectrics 300 (2004) 121. F. Shikanai, M. Komukae, Z. Czapla, T. Osaka, J. Phys. Soc. Jpn. 71 (2002) 498. H. Taniguchi, M. Machida, N. Koyano, J. Phys. Soc. Jpn. 72 (2003) 1111. E. Magome, M. Machida, Y. Tamura, M. Komukae, J. Phys. Soc. Jpn. 76 (2007) 124601. Ya. Shchur, A. Kityk, Phys. Status Solidi B 252 (2015) 476. I.R. Zachek, Ya. Shchur, R.R. Levitskii, O.B. Bilenka, Phys. B 452 (2014) 152. I.R. Zachek, Ya. Shchur, R.R. Levitskii, Phys. B 478 (2015) 113. J. Nayeem, T. Kikuta, N. Nakatani, F. Matsui, S.-N. Takeda, K. Hattori, H. Daimon, Ferroelectrics 332 (2006) 13. S. Dacko, Z. Czapla, J. Baran, M. Drozd, Phys. Lett. A 221 (1996) 217. M. Wiesner, Phys. Stat. Sol. (b) 238 (2003) 68. N. Yasuda, T. Sakurai, Z. Czapla, J. Phys.: Condens Matter 9 (1997) L347. J. Nayeem, H. Wakabayashi, T. Kikuta, T. Yamazaki, N. Nakatani, Ferroelectrics 269 (2002) 153. F. Shikanai, M. Yamasaki, M. Komukae, T. Osaka, J. Phys. Soc. Jpn. 72 (2003) 325.