Electrocaloric effect in KH$_2$PO$_4$ family crystals

The proton ordering model for the KH$_{2}$PO$_{4}$ type ferroelectrics is modified by taking into account the dependence of the effective dipole moments on the proton ordering parameter. Within the four-particle cluster approximation we calculate the crystal polarization and explore the electrocaloric effect. Smearing of the ferroelectric phase transition by a longitudinal electric field is described. A good agreement with experiment is obtained.


Introduction
The electrocaloric (EC) effect is the change of temperature of a dielectric at an adiabatic change of the applied electric field. Research in this field is driven by a quest for materials that can be used for efficient, environment-friendly, and compact (on-chip) solid-state cooling devices.
The current state of the art on the electrocaloric effect research for ferroelectrics is well summarized in [1,2]. At the moment, the largest effect is observed in perovskite ferroelectrics. Thus, in [3] in the PbZr 0.95 Ti 0.05 O 3 thin film with a thickness of 350 nm in a strong electric field (480 kV/cm) the obtained electrocaloric temperature change is ∆T = 12 K. Ab initio molecular dynamics calculations [4] predict ∆T ≈ 20 K in LiNbO 3 . In the hydrogen bonded ferroelectrics of the KH 2 PO 4 (KDP) type, the electrocaloric effect was studied for relatively low fields only. Thus, it has been obtained that ∆T ≈ 0.04 K at E ≈ 4 kV/cm [5], ∆T ≈ 1 K at E ≈ 12 kV/cm [6], and ∆T ≈ 0.25 K at T c and E ≈ 1.2 kV/cm [7].
Theoretical calculations of the electrocaloric effect in KDP have been made in [8] within the Slater model [9] and in the paraelectric phase only. It is also known that the Slater model gives incorrect results in the ferroelectric phase, and more complicated versions of the proton ordering model are required for an adequate description of these crystals. Thus, the effect of electric field on the physical characteristics of the KDP type crystals, such as polarization, dielectric permittivity, piezoelectric coefficients, elastic constants, has been described within the proton ordering model with the piezoelectric coupling to the shear strain ε 6 [10][11][12] and with proton tunneling [13] taken into account. However, these theories required, in particular, invoking two different values of the effective dipole moments for the paraelectric and ferroelectric phase [10,12]. This made impossible a correct description of the system behavior in the fields high enough to smear out the first order phase transition. There is an inner logical contradiction in the model: while no physical characteristic of a crystal should exhibit any discontinuity in the fields above the critical one, there is no smooth transition between the values of model parameters, rigidly set to be different for the two phases.
In the present paper we suggest a way to remove this contradiction. Assuming that the difference between the dipole moments is caused by non-zero values of the order parameter, we modify the proton ordering model accordingly. The field dependences of polarization, smearing of the first order phase transition, and the electrocaloric effect are described. direct lattice-mediated interactions [14,15]; σ q f is the operator of the z-component of a pseudospin, corresponding to the proton on the f -th hydrogen bond ( f = 1, 2, 3, 4) in the q-th cell. Its eigenvalues σ q f = ±1 are assigned to two equilibrium positions of a proton on this bond.
In (3),Ĥ sh is the Hamiltonian of short-range interactions between protons, which includes terms linear over the strain [11] Here, and ε, w, w 1 are the energies of proton configurations. The third term in (3) is a linear over the shear strain ε 6 field due to the piezoelectric coupling; ψ 6 is the deformational potential. The fourth term effectively describes the system interaction with the external electric field E 3 . Here, µ f is the effective dipole moment of the f -the hydrogen bond, and The fifth term in (3) is introduced in the present paper for the first time. It takes into account the assumed dependence of the effective dipole moment on the order parameter (pseudospin mean value) It is equivalent to a term proportional to P 3 3 E 3 in a phenomenological thermodynamic potential. Note that the terms like P 2 3 E 3 are not allowed because of the symmetry considerations, and we keep the Hamiltonian to be linear in the field E 3 .

43703-2
In view of the crystal structure of the KDP type ferroelectrics, the four-particle cluster approximation is most suitable for short-range interactions [15,16]. Long-range interactions and the termĤ E are taken into account in the mean field approximation. Thus, Combining the fourth term in (3) and the first term in (6), we obtain the following term in the Hamiltonian describes the jump of the dipole moment at the first order phase transition, its different values for the paraelectric and ferroelectric phase, and its smooth behavior in the fields above the critical one, when there is no jump of η. We can now use a single value of µ for both phases and remove the logical contradiction of the earlier theories, described in Introduction.
Proceeding with the standard calculations of the cluster approximation [10,12,16], we obtain the following expression for the proton ordering parameter The thermodynamic potential is then obtained in the following form Here, σ 6 is the formally introduced shear stress conjugate to the strain ε 6 . In numerical calculations we put σ 6 = 0. The condition of the thermodynamic potential minimum ∂G ∂ε 6 T,E 3 ,σ 6 = 0 yields an equation for the strain ε 6 In the same way, we derive the expressions for polarization P 3 and molar entropy of the proton subsys- Here, N A is the Avogadro number; R is the gas constant. The following notations are used: Expressions for dielectric susceptibilities, piezoelectric coefficients, and elastic constants derived [17] from equations (8), (9) are slightly different from the previous ones [10], where the dependence of the effective dipole moment on the order parameter was not taken into account. Numerical calculations, however, showed [17] that in zero electric field the difference is minor.
The molar specific heat of the subsystem described by the Hamiltonian (1) is Here, Notations introduced here are described in appendix. Then, the total specific heat is Here, ∆C σ is assumed to describe all the anomalies of the specific heat at the phase transition, whereas the regular background contribution to the specific heat, mostly from the lattice of heavy ions, is approximated by a linear temperature dependence As will be discussed later, this linear approximation agrees with the experimental data. Finally, the electrocaloric temperature change is calculated using the known formula where the pyroelectric coefficient is

Numerical calculations
To perform the numerical calculations we need to set the values of the following theory parameters: the Slater energies ε, w, w 1 ; the parameter of the long-range interactions ν c ; the effective dipole moment µ and the correction is due to proton ordering µ ′ ; the deformation potentials ψ 6 , δ s , δ a , δ 1 ;
They are chosen, obviously, by fitting the theoretical thermodynamic characteristics to the experimental data, as described in [12]. The obtained optimum sets of the model parameters are given in table 1.
To describe crystals with different deuteration levels, we use the mean crystal approximation, where the theory parameters are assumed to be linearly dependent on deuteron concentration (except for the parameter ν c , for which a small deviation from the linear dependence is assumed, as it is chosen from the condition that the calculated transition temperature coincides with the experimental one, which is also slightly non-linear). The dependence of the energy levels and interparticle interaction constants on deuteration is caused by the corresponding geometrical changes in the crystal structure with deuteration (elongation of the hydrogen bonds, changes in the distance between the equilibrium positions of H or D on the bonds, changes in the lattice constants, etc).   The primitive cell volume is taken to be v = 0.1946 · 10 −21 cm 3 for K(H 1−x D x ) 2 PO 4 and v = 0.202 · 10 −21 cm 3 for K(H 1−x D x ) 2 AsO 4 , irrespectively of the deuteration. The energy w 1 of proton configurations with four or zero protons near the given oxygen tetrahedron should be much higher than ε and w. Therefore, we take w 1 = ∞ (d = 0).
As we have already mentioned, when the dependence of the effective dipole moment on the order parameter is taken into account, the agreement between the theory and experiment for most of the calculated dielectric, piezoelectric, elastic characteristics, and specific heat of the studied crystals in the absence of an external electric field is neither improved nor worsened (see [17]). However, the present model allows us to describe more consistently the smearing of the first order phase in high electric fields.
The temperature dependence of the specific heat of KH 2 PO 4 and KD 2 PO 4 is shown in figure 1. The contribution ∆C σ is essential in the transition region and satisfactorily describes the experimental anomalies. As one can see, the total specific heat above T c can be well approximated by a linear temperature dependence, thus justifying the linear dependence of C regular , given by equation (14).
In figures 2 and 3 we plotted the temperature variation of polarization of K(H 1−x D x ) 2 PO 4 in different fields. The agreement with experiment is better at x = 0.89 (and 0.84, see [17]) than at x = 0. We believe this is due to proton tunnelling, essential in non-deuterated samples, which is not included in our model.
The field E 3 , which in these crystals is the field conjugate to the order parameter, induces non-zero polarization P 3 above the transition point. Polarization has a jump at T c , indicating the first order phase transition. With an increasing field, the polarization jump decreases, whereas the transition temperature  [19]; at x = 0.86 -△ [19]. Dashed lines 1' and 2': the theoretical results of [12]. and x = 0.89, respectively (c.f. 0.22 and 0.13 K cm/kV from our earlier calculations [10] and experimental 0.125 K cm/kV of [23] for x = 0.89). At some critical field E * , the jump vanishes, and the transition smears out. The calculated coordinates of the critical point are E * = 125 V/cm, T * c =122.244 K for x = 0 and 7.1 kV/cm, 212.55 K for x = 0.89, which agrees well with the experiment [22,23]. It should be noted that in our previous calculations [12] it was impossible to obtain a correct description of the polarization behavior in the fields above the critical one, due to the necessity of using two different values of the effective dipole moment µ in calculations.
The calculated electrocaloric changes of temperature ∆T of the K( crystals with the adiabatically applied electric field are shown in figures 4 and 5. The experimental data of [7] were obtained at T = 121 K, which was very close to the transition temperature of the sample used in the measurements. ♦ [20]. Symbols are experimental points; solid lines: the present theory; dashed lines: the theoretical results of [12]. the fields above 1 kV/cm are available. And, of course, possibilities for experimental measurements are limited by the dielectric strength of the samples. As one can see from the temperature dependence of ∆T (figure 6) for K(H 1−x D x ) 2 PO 4 crystals, the calculated electrocaloric temperature change is the largest at temperatures below T 0 c but close T c and can exceed 6 K; however, the fields required to reach ∆T that high are not accessible in reality, because most likely they exceed the dielectric strength of the crystals.

Conclusions
Taking into account the dependence of the effective dipole moment on the order parameter within the framework of the proton ordering model allows us to correctly describe the smearing of the ferroelectric phase transition in high electric fields as well as the electrocaloric effect in the KDP family crystals. The theory predicts the values of the electrocaloric temperature change of a few Kelvins in high fields.
Additional experimental measurements of ∆T in the fields above 2 kV/cm are necessary.