# Order-disorder, ferroelasticity and mobility of domain walls in multiferroic Cu-Cl boracite

C M Fernandez-Posada, Charlotte Cochard, J M Gregg, R W Whatmore, M

A Carpenter

► To cite this version:

C M Fernandez-Posada, Charlotte Cochard, J M Gregg, R W Whatmore, M A Carpenter. Order-disorder, ferroelasticity and mobility of domain walls in multiferroic Cu-Cl boracite. Journal of Physics: Condensed Matter, 2020, 33 (9), pp.095402. 10.1088/1361-648x/abcb0f. hal-05423582

HAL Id: hal-05423582

https://hal.science/hal-05423582v1

Submitted on 18 Dec 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
PAPER • OPEN ACCESS

# Order-disorder, ferroelasticity and mobility of domain walls in multiferroic Cu-Cl boracite

To cite this article: C M Fernandez-Posada et al 2021 *J. Phys.: Condens. Matter* **33** 095402

View the [article online](#) for updates and enhancements.

## You may also like

- *Strain behavior and lattice dynamics in Ni$_{2}$Mn$_{3}$In$_{7}$*
  C Salazar Mejía, A K Nayak, J A Schiemer et al.

- *Static and dynamic strain coupling behaviour of ferroic and multiferroic perovskites from resonant ultrasound spectroscopy*
  M A Carpenter

- *Magnetoelastic properties and behaviour of 4C pyrrhotite, Fe$_{7}$S$_{8}$, through the Besnus transition*
  C R S Haines, S E Dutton, M W R Volk et al.
# Order-disorder, ferroelasticity and mobility of domain walls in multiferroic Cu-Cl boracite

C M Fernandez-Posada<fnref n="1" />,*, C Cochard<fnref n="2" />, J M Gregg<fnref n="2" />,
R W Whatmore<fnref n="3" />,<fnref n="4" /> and M A Carpenter<fnref n="1" />

E-mail: cmf61@cam.ac.uk

Received 4 July 2020, revised 25 September 2020
Accepted for publication 17 November 2020
Published 11 December 2020

<figure><img src="image_2.png" /></figure>

CrossMark

## Abstract

Domain walls in Cu-Cl boracite develop as a consequence of an improper ferroelastic, improper ferroelectric transition, and have attracted close interest because some are conductive and all can be mechanically written and repositioned by application of an electric field. The phase transition and its associated dynamical properties have been analysed here from the perspective of strain and elasticity. Determination of spontaneous strains from published lattice parameter data has allowed the equilibrium long-range order parameter for $F43c \rightarrow Pca2_1$ to be modelled simply as being close to the order-disorder limit. High acoustic loss in the cubic phase, revealed by resonant ultrasound spectroscopy, is consistent with the presence of dynamical microdomains of the orthorhombic structure with relaxation times in the vicinity of $\sim 10^{-5}-10^{-6}$ s. Low acoustic loss in the stability field of the orthorhombic structure signifies, on the other hand, that ferroelastic twin walls which develop as a consequence of the order-disorder process are immobile on this time scale. A Debye loss peak accompanied by $\sim 1\%$ elastic stiffening at $\sim 40$ K is indicative of some freezing of defects which couple with strain or of some more intrinsic freezing process. The activation energy of $\ge \sim 0.01-0.02$ eV implies a mechanism which could involve strain relaxation clouds around local ferroelectric dipoles or freezing of polarons that determine the conductivity of twin walls.

Keywords: conductive domain walls, ferroelastic twin walls, phase transitions, boracite, multiferroic

(Some figures may appear in colour only in the online journal)

<figure><img src="image_3.png" /></figure>

Original content from this work may be used under the terms of the Creative Commons Attribution 4.0 licence. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI.

<footnote n="1">¹ Department of Earth Sciences, University of Cambridge, Downing Street, Cambridge CB2 3EQ, United Kingdom</footnote>

<footnote n="2">² Centre for Nanostructured Media, School of Mathematics and Physics, Queen's University Belfast,
Belfast BT7 1NN, United Kingdom</footnote>

<footnote n="3">³ Department of Chemistry, University College Cork, Cork, Ireland</footnote>

<footnote n="4">⁴ Department of Materials, Faculty of Engineering, Imperial College London, London, SW7 2AZ, United Kingdom</footnote>

<footnote>* Author to whom any correspondence should be addressed.</footnote>
# 1. Introduction

Boracites form an extensive family of compounds with general formula $M_3B_7O_{13}X$, where M is a divalent metal and X is generally a halogen [1, 2]. They possess a high-temperature cubic phase with space group $F\bar{4}3c$ and, depending on composition, undergo diverse phase transitions to produce structures with a range of different symmetries at room temperature, including tetragonal, trigonal, orthorhombic, or monoclinic [2, 3]. Essentially as a consequence of their diversity with respect to structure and symmetry, they have functional properties which are relevant for device applications. In addition to their piezoelectric properties [1], for example, their pyroelectric voltage figures of merit have made them interesting materials for use in infrared detectors [4, 5] and electro-optic applications [6]. Finding of the magnetoelectric effect and the presence of a toroidal moment in some boracites has also attracted attention [7] and effectively made them the first multiferroic compounds [8]. In addition to these impressive bulk properties, the most recent discovery of the existence of both electrically conducting and insulating charged domain walls (CDWs) makes them a serious prospect for use in new domain wall nanotechnologies [9].

Domain walls are regions of a crystal which facilitate the change of ferroic order parameters from one domain to another. They define steep gradients in strain, magnetic moment, ferroelectric dipole, or multiple combinations of these. As such, they can be considered as potential nanoscale functional materials in their own right, as they have local structure, symmetry and properties that necessarily differ from those of the domains which make up the bulk of the crystal [10–16]. Their unique properties open up possibilities for developing new devices based on tailored properties in the burgeoning fields of ‘domain boundary engineering’ [17, 18] and ‘domain wall nanoelectronics’ [19]. CDW’s provide one of the most interesting examples in this context [14, 16]. Electric dipoles organised head-to-head or tail-to-tail give a nominally charged wall [13, 16], and the existence of positive and negative bound charges at the walls creates an electrostatic potential that should be balanced by the re-distribution of carriers. Effectively, these walls carry bound charge while the rest of the material remains unaltered, resulting in quasi-two-dimensional charged sheets coupled to the neighbouring polarization. Ultrafast and ultralow energy consumption devices such as non-volatile memory could be based on the possibility of creating and erasing these walls, with control of data storage via the conductive channels [20].

Domain walls in Cu–Cl boracite, $Cu_3B_7O_{13}Cl$, are both ferroelastic and ferroelectric at room temperature and it has been shown that their configuration and functionality can be engineered by selective application of stress and electric fields [9]. Essential to such practical applications is an understanding of their mobility in response to external fields. For example, some devices will require that a particular configuration remains fixed for the lifetime of the device while high mobilities would be required if the objective is ultrafast switching. The primary objectives of the present study were to investigate,

specifically, the strain and elastic properties which accompany the cubic—orthorhombic transition in $Cu_3B_7O_{13}Cl$ and assess the dynamic properties of the domain walls which develop as a consequence of the symmetry change. Evolution of the symmetry-breaking shear strain reveals a temperature dependence for the driving order parameter which is consistent with the transition being predominantly order–disorder in character. Mechanical spectra obtained by resonant ultrasound spectroscopy (RUS) show that, in exact contrast with ferroelastic transitions in oxide perovskites, acoustic attenuation is high in the high symmetry (cubic) phase and low in the low symmetry (orthorhombic) phase. The attenuation at high temperatures is attributed to the existence of locally ordered microdomains. Low attenuation below the transition point shows that the domain walls themselves are immobile under the low stress and high frequency conditions that apply in RUS experiments. By analogy with the properties of domain walls in perovskites, the mobility might be substantially enhanced at low temperatures in boracites which are proximal to additional structural, magnetic or electronic instabilities. The low activation energy of an additional Debye-like loss peak at ~40 K reveals a further freezing process.

# 2. Order/disorder character and strain coupling

The approach of Dvorak et al [21, 22] in using symmetry and coupling of physical properties with the order parameter for characterising phase transitions in boracites is used here as a starting point. As set out in the appendix, the active representation for the $F\bar{4}3c-Pca2_1$ transition at ~365 K in $Cu_3B_7O_{13}Cl$ is $X_5$ and the order parameter has direction P10 [21–23]. Equation (A1) is a Landau expansion in one order parameter, $q$, including coupling with polarisation, $P$, and symmetry breaking shear strains, $e_t$ and $e_6$, of the form $\lambda Pq^2$, $\lambda e_t q^2$ and $\lambda e_6 q^2$. It follows that $P, e_t, e_6$ are expected to vary with $q^2$ and, hence, that the transition is improper ferroelectric, improper ferroelastic. The structural changes involved are illustrated in figure 1(a). A dipole oriented in the direction [001]$_{pc}$ is created and there are six equivalent twin variants [figure 1(b)]. (Subscript pc refers to the cubic parent structure). Details of these twins and their influence on Raman spectra have been set out in Iliev et al [24].

The evolution of the order parameter can be determined through analysis of the spontaneous shear strain and comparison with data for polarisation and birefringence. Spontaneous strains calculated from lattice parameters reported by Uesu et al [25] in the range 200–450 K are set out in the appendix with a choice of reference axes such that the symmetry-breaking shear strains are $e_t = \frac{1}{\sqrt{3}}(2e_3 - e_1 - e_2)$ and $e_o = (e_1 - e_2)$, rather than $e_t$ and $e_6$. The dominant shear strain is $e_o$ and figure 2 shows that this evolves with temperature in the same manner as the measured values of optical birefringence in the (110) plane, $\Delta n$ [26], and polarisation [5]. Close correlation between these properties confirms that they provide an internally consistent description of the evolution of $q$ as $e_o \propto P \propto \Delta n_{(110)} \propto q^2$. Also shown is the excess entropy calculated from heat capacity data [27] using the approach
<figure><img src="image_4.png" /><figcaption>**Figure 1.** (a) Schematic representation of the transformation of Cu₃B₇O₁₃Cl from the high-temperature cubic F̄43c structure to the orthorhombic Pca2₁ structure. The projections show the similarity of both crystal structures. Cl atoms appear to be hidden in the cubic F̄43c projection. (b) Three of the six twin variants in the orthorhombic structure. The remaining three are created by exchanging the [100]₀ and [010]₀ axes.</figcaption></figure>

set out, for example, in Carpenter *et al* [28]. The transition is weakly first order at $T_{tr} = 365$ K.

Evidence from inelastic neutron scattering studies is that cubic Cu–Cl boracite does not have a soft mode at the X point of the Brillouin zone and that the high temperature structure can be represented as containing long-lived microdomains of the orthorhombic structure [29]. It is therefore possible that the mechanism is of order/disorder rather than displacive character. Issues of displacive and order/disorder contributions (eg in references [30, 31]) can be tested at the simplest phenomenological level by comparison of Landau and Bragg–Williams solutions for the temperature dependence of the order parameter, *q*. As set out in the appendix, neither of the standard Landau solution for a weakly first order transition, representing a displacive limit, or the standard Bragg–Williams solution [32], representing an order/disorder limit, reproduces the temperature dependence of *q*². It turns out, however, that only a small modification of the Bragg–Williams equation is needed to provide a satisfactory representation of the experimental data, apart from in the temperature interval immediately below the transition point where dynamical effects may still be important. The fit shown in figure 3 is the solution to the excess free energy for a 2-site order/disorder model in which a fourth order term in *q* is included, following Malcherek *et al* [33]:

$$G = \frac{1}{2}Aq^2 + \frac{1}{4}Bq^4 - TS_{\text{config}}. \quad (1)$$

The configurational entropy is:

$$S_{\text{config}} = -R[(1+q)\ln(1+q) + (1-q)\ln(1-q)]. \quad (2)$$

Values of the fit parameters, $A = -5.75$ kJ mol⁻¹ and $B = -4.06$ kJ mol⁻¹, give a first-order transition with a discontinuity between $q = 0$ and $q = 0.69$ at the transition point, $T_{tr} = 374$ K. Failure to match the exact transition point and evolution of the order parameter in the ~40 K interval below it is not surprizing given that this is a static model. A more complete description would need to include increasing dynamical contributions as the transition point is approached.

Experimental values of the excess entropy obtained by integration of the excess heat capacity reach ~9 J (mol⁻¹ K⁻¹) at 200 K [figure 2(c)] and are consistent with a total configurational entropy change for $q = 0 \rightarrow q = 1$ of –R2 ln 2 = –11.5 J (mol⁻¹ K⁻¹) expected on the basis of equation (2). Equation (2) also gives an approximately linear dependence of ΔS on q² up to q ~0.7–0.8 for an order/disorder process [32] so that all the correlations between measured properties in figure 2 are consistent with a simple description of the transition as being close to the order/disorder limit.

The magnitudes of strains coupled with the order parameter vary up to a few % [appendix figure A1(b)] and are thus comparable with the strains which accompany octahedral tilting transitions in perovskites [34, 35].

### 3. Sample description and experimental details

The small polished crystal used for determination of strain-related dynamic behaviour associated with ferroelasticity in Cu₃B₇O₁₃Cl was prepared by the phase transport technique [36] and came from the same single crystal as used in the work of McQuaid *et al* [9]. Its shape was close to that of a rectangular
<figure><img src="image_5.png" /><figcaption>**Figure 2.** Orthorhombic shear strain, $e_o (\propto q^2)$, plotted against temperature showing scaling with (a) spontaneous polarization, *P* (data from [5]), (b) birefringence in (110), $\Delta n(110)$ (data from [26]), and (c) excess entropy, $\Delta S$ (from integration of excess heat capacity data given by Delfino *et al* [27]). Note: experimental uncertainties were not given in the original papers. Some indication of relative uncertainties in each dataset is provided by the magnitude of scatter between data point.</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>**Figure 3.** Evolution of the order parameter, *q*, from the solution to equation (1), as fit to the evolution of the square root of the orthorhombic strain. The solution shown has a discontinuity between *q* = 0 and *q* = 0.69 at $T_{tr}$ = 374 K. The fit is less good in the immediate temperature interval below the transition point, where dynamical effects are likely to be important.</figcaption></figure>

parallelepiped, with dimensions ~1 × 2 × 0.29 mm³ and mass 0.0024 g. The largest pair of faces were aligned parallel to (100)$_{pc}$. These faces were polished so that configurations of twin domains could be observed between crossed polars in a transmitted light microscope.

The RUS technique for measuring elastic and anelastic properties of mm-sized samples has been described in detail by Migliori and Sarrao [37]. In essence, a piezoelectric transducer is used to excite natural resonances of the sample which are then detected by a second transducer. Among other applications, RUS has been used to determine the form and strength of coupling of the order parameter with strain at ferroic phase transitions and to investigate the dynamics of ferroelastic twin wall motion [38, 39]. Three different instruments were used in

the present study. High-temperature spectra were collected at temperatures between room temperature and ~450 K with the crystal sitting between the tips of alumina buffer rods inserted into a horizontal furnace [39]. The driving and detecting transducers were attached to the ends of the buffer rods outside the furnace, and a maximum voltage of 1 V was applied to the driving transducer. Note that the upper temperature limit of ~450 K was chosen to ensure that there was no possibility of thermal decomposition of the sample.

Low-temperature RUS measurements were made in a helium-flow cryostat [40] and in an Oxford instruments Teslatron cryostat [41]. The same design of head was used in the two low temperature instruments, with the sample held directly between the two transducers. Maximum voltages applied to the driving transducer were 10 V and 2 V, respectively.

Automated heating/cooling cycles included a dwell time of 15 min before the start of data collection at each fixed temperature, to allow for thermal equilibration. Each spectrum contained 65 000 data points in the frequency range 0.1–1.2 MHz. Individual resonance peaks were fit with an asymmetric Lorentzian function to determine their peak frequency, *f*, and width at half maximum height, $\Delta f$, using the software package IGOR (wavemetrics). Values of $f^2$ for different peaks scale with different combinations of predominantly shear elastic constants, as illustrated for example by the evolution of the shear modulus of polycrystalline samples [35, 42]. Acoustic loss is expressed in terms of the inverse mechanical quality factor, $Q^{-1}$, and in RUS experiments is generally specified as being equal to $\Delta f/f$. The loss mechanisms are likely to be dominated by anelastic relaxation of defects in the crystal, such as ferroelastic twin walls, in response to an applied dynamic shear stress for which the theory is set out in full by Nowick and Berry [43].

A quantum design PPMS dynacool instrument with 9T superconducting magnet and the ACMS II option was used to carry out magnetic measurements on the same crystal as used for RUS. DC measurements were made in two different configurations: (1) with the (100)$_{pc}$ faces of the sample perpendicular to the magnetic field, and (2) with the sample mounted on a
<figure><img src="image_7.png" /><figcaption>**Figure 4.** (a) Stack of segments of RUS spectra collected during heating and cooling in steps of ~2 K through the temperature interval 290–454 K. Resonance peaks from the sample show marked softening with increasing temperature and are absent from spectra collected at ~366 K and above, which is immediately above the cubic–orthorhombic transition point. Resonance peaks which do not display any change in frequency as a function of temperature, marked with a star, are from the buffer rods. Spectra collected at 365 K are highlighted in green (cooling) and black (heating). (b) Variation with temperature of $f^2$ (circles) and $Q^{-1}$ (crosses) from fitting of a resonance peak which had frequency ~730 kHz at room temperature. (c) and (d) Ferroelastic twins in the Cu–Cl crystal, as observed at room temperature between crossed polars in an optical microscope before (c) and after (d) heating through the cubic–orthorhombic transition. Scale bars = 750 µm.</figcaption></figure>

quartz holder so that the (100)$_{pc}$ faces were parallel to the field. Subscript pc refers to the cubic reference state. Data were collected during heating through the temperature interval 2–50 K in a measuring field of 1000 Oe following cooling in zero field (ZFC) or following cooling in a 1000 Oe field (FC). AC measurements through the temperature interval 2–60 K were also made with the sample held with its (100)$_{pc}$ faces perpendicular to the field using a DC field of 20 Oe and an AC amplitude of 5 Oe at frequencies of 100 Hz, 1 kHz and 10 kHz.

# 4. Results

## 4.1. Elasticity and anelasticity at high temperatures

Figure 4(a) shows segments of RUS spectra collected during heating and cooling through the interval 290–454 K. The y-axis represents amplitude but, in order to display the temperature dependence more clearly, the spectra have been stacked in proportion to the temperature at which they were collected and the axis relabelled as temperature. Peaks showing marked temperature-dependence are from the sample while peaks which are independent of temperature are from the buffer rods. In spectra collected within the stability field of the orthorhombic structure at ~363 K and below, the former are either clearly visible between the latter or are evident from enhanced amplitude where the two overlap. The cubic–orthorhombic transition expected at ~365 K is marked

by the complete disappearance of resonance peaks from the sample, indicating a high degree of acoustic attenuation in the cubic phase. This is particularly striking as the pattern is the opposite of what would usually be expected. Previous experience of many ferroelastic phase transitions is that the acoustic loss is larger in the low-temperature phase because the loss mechanism typically involves the mobility under external stress of ferroelastic twin walls [38]. Meyer *et al* [29] found evidence from neutron diffraction data for the presence of microdomains of orthorhombic structure in cubic Cu–Cl boracite and the loss process in the high-temperature phase is most likely due to dynamical disordering of the microdomains on the resonance time scale of ~$10^{-5}$–$10^{-6}$ s.

Variations of $f^2$ and $Q^{-1}$ from fitting of a resonance peak with frequency near 730 kHz at room temperature [figure 4(b)] reveal in more detail how the elastic properties evolve below $T_{tr}$. Values of $Q^{-1}$ immediately below the cubic–orthorhombic transition are ~0.002–0.003 and there is perhaps a slight tendency for these to increase with falling temperature below ~310 K. It is clear, however, that acoustic attenuation associated with any mobility of the twin walls is low. There are small differences in the values of $f^2$ between heating and cooling which are most likely to have been due to a change in the configuration of ferroelastic twins in the crystal. Such a change is confirmed by differences between ‘before’ and ‘after’ images of the crystal obtained between crossed polars in a transmitted light microscope [figures 4(c) and (d)].
<figure><img src="image_8.png" /><figcaption>**Figure 5.** (a) Stack of segments of RUS spectra showing a resonance peak with frequency near 660 kHz. The spectra were obtained during cooling from 297 K. (Teslatron cryostat, 2 V applied to driving transducer) (b) $f^2$ (circles) and $Q^{-1}$ (crosses) variations as a function of temperature for the same resonance peak in the full cooling/heating cycle between 299 and 4 K. The solid black curve is a fit of equation [3] to $Q^{-1}$ data in the interval 10–75 K. Fit parameters: $E_a/r_2(\beta) = 2.2 \pm 0.3$ kJ mol$^{-1}$, $T_m = 45.5 \pm 0.9$ K, $Q_m = 0.0008 \pm 0.0001$. The linear baseline is shown as a black line (c) $f^2$ (circles) and $Q^{-1}$ (crosses) variations as a function of temperature for a repeat cooling/heating cycle. Fit parameters for the Debye peak: $E_a/r_2(\beta) = 0.91 \pm 0.06$ kJ mol$^{-1}$, $T_m = 35.4 \pm 0.7$ K, $Q_m = 0.0014 \pm 0.0001$ (d) $f^2$ and $Q^{-1}$ variations as a function of temperature for a cooling/heating cycle in the He-flow cryostat with 10 V as the maximum voltage applied to the driving transducer. $Q^{-1}$ values increased by a factor of ~2 due to the increase in applied voltage.</figcaption></figure>

## 4.2. Elasticity and anelasticity at low temperatures

Figure 5(a) shows a single resonance peak with frequency near 660 kHz in segments of RUS spectra stacked up the y-axis according to the temperature at which they were collected during heating from 6 to 297 K in steps of ~5 K in the Teslatron cryostat. Figure 5(b) shows the evolution of $f^2$ and $Q^{-1}$ from fitting of the peak for the full cooling/heating cycle. There is a small anomaly in the temperature dependence at ~150 K which has not been accounted for, but the clearest feature is a small degree of stiffening at ~50 K followed by softening below ~25 K. The Debye-like peak in $Q^{-1}$ is centred on ~40 K, which coincides with an increase in $f^2$ by ~1%. To confirm this low temperature pattern, measurements were repeated to include smaller temperature steps below 100 K [figure 5(c)]. There is a trend of increasing $Q^{-1}$ with increasing temperature above ~100 K but the anomaly at ~150 K does not appear to have been reproducible. There are also small differences between $f^2$ values obtained during cooling and heating.

A thermally activated acoustic loss process followed as a function of temperature at more or less constant frequency can be described by [44–47]

$$Q^{-1}(T) = Q_m^{-1} \left[ \cosh \left\{ \frac{E_a}{Rr_2(\beta)} \left( \frac{1}{T} - \frac{1}{T_m} \right) \right\} \right]^{-1} . \quad (3)$$

The maximum value of $Q^{-1}$, $Q_m^{-1}$, occurs at temperature $T_m$ and depends on the difference in elastic constants between the relaxed and unrelaxed states. The width of the peak depends on both the activation energy, $E_a$, and the parameter $r_2(\beta)$ which relates to the spread of relaxation times. If the dissipation process involves a single relaxation time, the value of $\beta$ is 0 and the value of $r_2(\beta)$ is 1. R is the gas constant. The fit of equation (3) shown in figure 5(b) includes a linear baseline and has $E_a/r_2(\beta) = 2.2 \pm 0.3$ kJ mol$^{-1}$ (~0.02 eV), $T_m = 45.5 \pm 0.9$ K, $Q_m^{-1} = 0.0008 \pm 0.0001$. The fitting gave an adjusted $R^2$ value of 0.915 and $\chi^2 = 7.38 \times 10^{-9}$. Larger values of $r_2(\beta)$ would apply if there was a spread of relaxation times so ~0.02 eV is a lower limit for the true activation energy. The fit of equation (3) shown in figure 5(c) has $E_a/r_2(\beta) = 0.91 \pm 0.06$ kJ mol$^{-1}$ (~0.01 eV), $T_m = 35.4 \pm 0.7$ K, $Q_m^{-1} = 0.0014 \pm 0.0001$ ($R^2 = 0.807$, $\chi^2 = 4.43 \times 10^{-8}$). Assuming Arrhenius behaviour, $f = f_o \exp(-E_a/RT)$, allows values of the constant, $f_o$, to be determined as $1.46 \times 10^7$ Hz for $E_a = 0.91$ kJ mol$^{-1}$ and $2.23 \times 10^8$ Hz for $E_a = 2.2$ kJ mol$^{-1}$. For comparison with magnetic data below, these would give $T_m$ values of 11.4 K and 21.4 K at 1 kHz, or 15.9 and 26.4 K at 10 kHz.
<figure><img src="image_9.png" /><figcaption>**Figure 6.** Temperature dependence of the field cooled DC magnetization, *M*, in the temperature range 2–50 K, measured in a 1000 Oe field (a) with the field parallel and (b) perpendicular to the (100)\*<sub>pc</sub> faces. ZFC stands for zero-field cooled and FC for field cooled.</figcaption></figure>

A loss peak in acoustic data implies freezing of some dynamical process in the sample which is coupled with strain. In relation to the possible contributions of defects, an activation energy of ~0.01 eV and elastic stiffening by only ~1% rules out a mechanism involving pinning of ferroelastic twin walls by oxygen vacancies. The activation energy for this is on the order of ~0.5–1 eV in perovskites, for example [38]. A more likely loss process could involve motion under external stress of small local strain fields around magnetic or ferroelectric defects. For example, activation energies of ~0.05 eV fall within a range that might be associated with freezing of polarons [41, 45]. One interesting possibility, in this context, is that the freezing process is of polarons responsible for conductivity at the ferroelastic/ferroelectric twin walls. More collective behaviour, such as the formation of polar nanoregions in relaxor ferroelectrics, can also have activation energies of about the same order of magnitude [42]. The most obvious correlation with other measured properties reported in the literature for Cu₃B₇O₁₃Cl is a broad hump in magnetic susceptibility, $\chi(T)$, which has its maximum value at ~20 K [3]. This is ahead of the magnetic transitions reported to occur at ~8 K and ~9 K [3].

RUS measurements repeated in the He-flow cryostat using 10 V as the maximum voltage applied to the driving transducer revealed an additional feature of the loss mechanism. Values of $Q^{-1}$ extracted from fitting of the resonance peak

<figure><img src="image_10.png" /><figcaption>**Figure 7.** (a) Real ($\chi'$) and (b) imaginary ($\chi''$) parts of ac magnetic susceptibility from 2 to 60 K, with an external dc magnetic field of 20 Oe and an ac field with amplitude 5 Oe at 1 kHz and 10 kHz from a sample with mass 0.0024 g.</figcaption></figure>

near 660 kHz increased by a factor of ~2 and values of $f^2$ by ~1%–2% across the full temperature range [figure 5(d)], in comparison with the earlier measurements at 2 V. The loss peak at ~40 K became less clearly defined and the peak in $f^2$ at ~20 K became smoothed out. A repeat set of measurements and fitting of the peak with frequency near 760 kHz revealed that the low-temperature anomalies are not suppressed for all resonance modes, but increasing the voltage in this way results in increases in amplitude of resonance peaks. The increase in stress required to increase the strain amplitude is evidently sufficient to reduce pinning/freezing of the defect responsible for the Debye loss peak.

### 4.3. DC and AC magnetism

Magnetic measurements were carried out in order to test whether the Debye loss peak that appears in the RUS measurements correlates with any magnetic process. Figure 6 contains the DC magnetization, *M*, curves for the sample when the (100)\*<sub>pc</sub> faces were perpendicular [figure 6(a)] and parallel [figure 6(b)] to the applied field. Two small anomalies evident in data for (100)\*<sub>pc</sub> perpendicular to the applied field correspond closely with previously reported magnetic transitions which were presumed to be from $Pca2_11'(mm21')$ to an mm2 structure ($T_{N2}$ = 8.9 K here) and then from mm2 to $Pc'a'2_1(m';m'2)$ ($T_{N1}$ = 7.7 K here) [3]. The additional broad peak near 20 K was previously suggested to be characteristic of low-dimensional quantum magnets and/or frustrated spin systems [3]. When the (100)\*<sub>pc</sub> face was parallel to the field, only
the transition at 8.9 K and the hump near 20 K were visible
[figure 6(b)].

AC measurements were carried out to test whether there
is a dynamic magnetic component to the broad peak in the
DC data. Figures 7(a) and (b) show respectively the real (χ')
and imaginary (χ'') parts of AC susceptibility for the set up in
which the applied field was parallel to the (100)_{pc} face of the
sample. They were collected in a heating sequence, following
cooling in zero-field. Only the results for 1 kHz and 10 kHz
are shown because the data collected at 100 Hz were too weak
to be informative. The broad hump in χ' correlates with the
broad peak in DC magnetisation but is not overly dependent
on frequency. χ'' perhaps has some small anomaly near ~30 K
but also does not display any detectable difference between 1
and 10 kHz. The signal obtained with the field perpendicular
to the (100)_{pc} face was too weak to be informative at all three
measuring frequencies.

On the basis of the estimated kinetic parameters given
above for the acoustic loss peak measured at ~660 kHz, a sim-
ilar dynamic effect would be expected at temperatures of ~15
or ~25 K at measuring frequencies of 1 or 10 kHz. This falls
in the temperature interval of the broad peak in DC and AC
magnetic susceptibility. However, there is no evidence of dis-
persion with respect to frequency in the AC data so there is
not an overt correspondence between the mechanism respon-
sible for the acoustic loss peak and the origin of the magnetic
anomaly.

5. Discussion and conclusions: acoustic
properties of domain walls in a multiferroic
order/disorder system

Analysis of the cubic–orthorhombic transition in Cu–Cl
boracite from the perspectives of strain and elasticity has
revealed that it can be understood in terms of a driving
mechanism which is order/disorder rather than displacive in
character. The existence of dynamical microdomains of the
orthorhombic structure in the stability of the cubic structure
appears to be confirmed by complete attenuation of acoustic
resonances on timescales of ~10<sup>-5</sup>–10<sup>-6</sup> s. The attenuation
persisted up to the highest temperature at which data were col-
lected, ~450 K, which is substantially above the temperature
of T<sub>tr</sub> + 10 K where more direct evidence of microdomains
was provided by inelastic neutron scattering [29]. An increase
in acoustic loss above the transition point also mirrors the
behaviour of some metal-organic frameworks, in which the
transition to a higher symmetry state allows dynamical dis-
ordering of the organic molecule [48–50]. This pattern of
attenuation is quite different from that typically observed at
displacive phase transitions, such as due to octahedral tilt-
ing in perovskites for example, where attenuation in the high
symmetry parent crystal is low [38].

The magnitudes of symmetry breaking shear strains iden-
tified here for Cu–Cl boracite are comparable with those
which accompany octahedral tilting transitions in oxide per-
ovskites such as (Ca,Sr)TiO₃, PrAlO₃, LaAlO₃ and SrZrO₃
[34, 35, 51–53]. All other things being equal, the mobility
of domain walls might be expected also to be comparable,

therefore. There is no evidence, as yet, that an order/disorder
mechanism for the transition necessarily results in slower
domain wall motion and/or stronger pinning by local strain
fields. Attenuation sufficient to cause the disappearance of res-
onance peaks in RUS spectra (superattenuation) in ferroelastic
phases of (Pr,La)O₃ and (Ca,Sr)TiO₃ occurs in the tempera-
ture interval ~900–500 K [54–56], below which the domain
walls become pinned by oxygen vacancies. Essentially the
same tilting transition occurs in EuTiO₃ at 284 K, which is
well below the expected freezing/pinning temperature of the
domain walls, and there is only a relatively small peak in Q⁻¹
immediately below the transition point followed by a steady
decrease with falling temperature [57, 58]. The simplest expla-
nation for the drop in Q⁻¹ values below Tₜr observed here
is then that the pinning/freezing temperature for the domain
walls at measuring frequencies near 1 MHz is above 365 K.
It is likely that the walls can be unpinned at low temperatures
by increasing the stress above some critical value, as has been
observed in LaAlO₃ [59].

An important factor which affects domain wall dynamics is
the width of the walls, since thin walls are more likely to be
pinned by point defects than thick walls. According to Catalan
et al [60], domain walls which are both ferroelectric and fer-
roelastic are expected to be thinner than pure ferroelastic walls
but thicker than purely ferroelectric walls. On the other hand,
Jia et al [61] found that charged twin walls in PbZr<sub>0.2</sub>Ti<sub>0.8</sub>O<sub>3</sub>
thin films are ~10 times thicker than uncharged walls. A
change in the internal structure of twin walls in SrTiO<sub>3</sub>, result-
ing in the development of polarity due to local ferroelectric
displacements, probably accounts for their high mobility down
to at least ~5 K [62]. Ferroelastic twin walls due to octahedral
tilting in PrAlO<sub>3</sub> also cause superattenuation of acoustic reso-
nances down to at least 10 K when a second (electronic) order
parameter develops [63]. It remains to be seen whether high
domain wall mobilities might be induced in boracites which
undergo additional structural or electronic instabilities below
room temperature but such a change would be a signal of local
structure and/or thickness of the walls.

Fitting of the Debye-like peak in Q⁻¹ observed at ~40 K
[figures 5(b) and (c)] gave a low activation energy, consis-
tent with a loss mechanism related to freezing of polar mag-
netic/ferroelectric defects or of polarons related to the dis-
appearance of conductivity at domain walls. Increasing the
driving voltage to increase the amplitude of the intrinsic res-
onances of the sample was sufficient to substantially reduce
the attenuation, consistent with relatively weak pinning. Addi-
tional measurements of the bulk dielectric constant and elec-
trical conductivity of domain walls at these temperatures are
needed to narrow down on possibilities for the loss mechanism
but, on the basis of the AC magnetic data at least, it does not
appear to have a strong magnetic component. There is, again,
a close analogy with oxide perovskites in that the loss peak
is similar to a peak in Q⁻¹ seen at ~50 K in RUS data from
LaAlO₃ [64].

Finally, there is known to be significant magnetoelectric
coupling in Cu–Cl boracite, with the ferromagnetic moment
and spontaneous polarisation both aligned parallel to [001]pc
of the orthorhombic structure below ~8–9 K [3, 65, 66]. As
a consequence, the ferroelastic/ferroelectric twin walls must also contain gradients of magnetic moment and it is almost inevitable that the temperature of the magnetic ordering transition would be slightly different at the domain walls in comparison with inside the domains. It is worth noting that heat capacity data of Schnelle and Schmid [3] and magnetic data in figure 6(a) show evidence of two closely spaced transitions. One possible interpretation is that these correspond to separate transitions within the domain walls and within the domains.

In conclusion, although the twin walls of Cu–Cl boracite are multiferroic, in the sense that they must inevitably have steep gradients of electric dipole moment and ferroelastic shear strain, and the transition mechanism is order/disorder rather than displacive, there is as yet no evidence that their response to applied stress is fundamentally different from that of ferroelastic domain walls in many perovskites. Conversion of theoretical possibilities for using the conductive properties of the domain walls in practical devices will require that they can be manipulated into some preferred configuration using external fields and then held in place for the lifetime of the device. It should be possible to engineer the balance between slow and fast responses in exactly the same manner as might be used for ferroelastic domain walls due to tilting transitions in oxide perovskites. The key pinning mechanisms are likely to involve strain coupling with point defects in the underlying crystal structure and it should be possible to increase their mobility by choosing a boracite phase which is close to some additional structural, magnetic or electronic instability. In this context, choice of the magnetic cation should provide particular interesting possibilities for inducing phenomenologically rich behaviour at the domain walls.

## Acknowledgments

This work was funded by EPSRC Grant No. EP/P024904/1 and US-Ireland R&D Partnership Programme (Grant No. USI 120). RUS facilities were established through grants from the Natural Environment Research Council (Grants Nos. NE/B505738/1 and NE/F017081/1) and the Engineering and Physical Sciences Research Council (Grant No. EP/I036079/1) to MAC. We thank Cheng Liu for the magnetic measurements carried out at Maxwell Centre, Cavendish Laboratory.

## Appendix. Strain analysis

A Landau expansion obtained using the software package ISOTROPY [67] for the improper ferroelastic/ferroelectric transition $F\bar{4}3c \rightarrow Pca2_1$, driven by an order parameter, $q$, which has the symmetry of irreducible representation $X_5$ (direction P10: $(0,0,0,0,0,a,-a)$) coupled with polarisation, $P$, ($\Gamma_4$, direction P1(3): $(0,0,b)$, polarisation parallel to $[001]_{pc}$), and strain, $e$, is:

$$G = \frac{1}{2}a(T - T_c)q^2 + \frac{1}{4}bq^4 + \frac{1}{6}cq^6 + \lambda_1 e_v q^2 + \lambda_2 e_t q^2 + \lambda_3 e_6 q^2 + \lambda_4 e_0^2 q^2 + \lambda_5 (e_4^2 + e_5^2) q^2 + \lambda_6 e_4 e_5 q^2 + \lambda_7 e_6 P$$

$$+ \lambda_8 e_v P^2 + \lambda_9 e_t P^2 + \lambda_{10} e_4 e_5 P + \lambda_{11} e_0^2 P^2 + \lambda_{12} q^2 P + \frac{1}{4} (C_{11}^0 - C_{12}^0) (e_0^2 + e_t^2) + \frac{1}{6} (C_{11}^0 + 2C_{12}^0) e_v^2 + \frac{1}{2} C_{44}^0 (e_4^2 + e_5^2 + e_6^2). \quad (A1)$$

$G$ is the excess free energy due to the transition, $a, b, c$ are standard Landau coefficients, $T_c$ is the critical temperature, $e_v$ is volume strain, $e_0$ and $e_t$ are orthorhombic and tetragonal shear strains defined with respect to reference axes which are parallel to the crystallographic axes of the cubic structure, $e_4, e_5, e_6$ are shear strains, $\lambda_i$ are coupling coefficients and $C_{ik}^0$ are elastic constants of the cubic reference structure. The standard solution for a first order transition, such as the $\alpha-\beta$ transition in quartz [28] is:

$$q^2 = \frac{2}{3} q_0^2 \left\{ 1 + \left[ 1 - \frac{3}{4} \left( \frac{T - T_c}{T_{tr} - T_c} \right) \right]^{1/2} \right\}, \quad (A2)$$

where $q_0$ is magnitude of the discontinuity in $q$ at the transition temperature, $T_{tr}$.

For convenience, reference axes $X$ and $Y$ used to calculate the spontaneous strains have been rotated through $45^\circ$ in comparison with the reference system used in ISOTROPY. This means that $e_6$ and $e_0$ are exchanged, leading to the expected relationships between $q$ and the symmetry-adapted strains as $e_v \propto e_0 \propto e_t \propto q^2$, with $e_4 = e_5 = e_6 = 0$. Values of spontaneous strains $e_1, e_2, e_3, e_0, e_t$ have been calculated in the usual way [68] using:

$$e_1 = \frac{a - a_0}{a_0} \quad (A3)$$

$$e_2 = \frac{b - a_0}{a_0} \quad (A4)$$

$$e_3 = \frac{c - c_0}{c_0} \quad (A5)$$

$$e_v = e_1 + e_2 + e_3 \quad (A6)$$

$$e_0 = e_1 - e_2 \quad (A7)$$

$$e_t = \frac{1}{\sqrt{3}} (2e_3 - e_1 - e_2). \quad (A8)$$

$a, b$ and $c$ are the lattice parameters of the orthorhombic structure given by Uesu et al [25] and $a_0$ and $c_0$ are reference lattice parameters of the cubic phase extrapolated into the stability field of the orthorhombic phase. The cubic structure should have $a_0 = c_0/\sqrt{2}$, but the way in which the original x-ray data were collected resulted in a slight discrepancy from this relationship. The discrepancy does not affect the strain analysis to any significant extent.

Figure A1(a) contains lattice parameters of Cu–Cl boracite reproduced from Uesu et al [25] and includes straight line fits to data for the cubic phase, with extrapolation to give $a_0$ and $c_0$. Values of the linear strains and the volume strain are given in figure A1(b) and the evolution of the symmetry-breaking shear strains, $e_0$ and $e_t$, are given in figure A1(c). The dominant shear
<figure><img src="image_11.png" /><figcaption>**Figure A1.** (a) Lattice parameter data for Cu–Cl boracite reproduced from Uesu *et al* [25], as specified with respect to crystallographic axes of the orthorhombic structure. Linear fits to data for the cubic phase are $a_o = 8.45 + 5.85 \times 10^{-4}T$, $c_o = 11.95 + 6.87 \times 10^{-5}T$. (b) Strains $e_1$, $e_2$, $e_3$ and $e_v$ as a function of temperature, calculated from equations (A3)–(A6). (c) Orthorhombic strain ($e_o$) and tetragonal strain ($e_t$) as function of temperature, given by equations (A7) and (A8). Note: experimental uncertainties were not given in the original work.</figcaption></figure>

strain is $e_o$, and all the strain parameter show a discontinuity at $T_{tr} \approx 365$ K.

The variation of $e_o$ can be used to test models for the temperature dependence of $q^2$, on the basis that the expected relationship is $e_o \propto q^2$. Equation (A2) does not provide a good description of this evolution over the full temperature range and fitting to two separate ranges (figure A2) is also inadequate. Thus the Landau solution for a weakly first order transition describes the transition only on a qualitative basis.

The standard solution for the Bragg–Williams model of a thermodynamically continuous, pure order/disorder transition occurring at $T_c$, as reproduced in [32], is

<figure><img src="image_12.png" /><figcaption>**Figure A2.** Evolution of the orthorhombic strain ($e_o \propto q^2$) with temperature, together with fits of equation (A2) in two different ranges.</figcaption></figure>

<figure><img src="image_13.png" /><figcaption>**Figure A3.** Evolution of the order parameter, $q$, obtained from the Bragg–Williams solution, equation (A9), in comparison with the evolution of the orthorhombic strain, $e_o$, taking into account $e_o \propto q^2$.</figcaption></figure>

$$q = \tanh \left( \frac{T_c}{T} q \right). \qquad (A9)$$

This has been compared to the data for $\sqrt{e_o}$ in figure A3. As with the Landau solution, the fit is only qualitative at best.

## ORCID iDs

C M Fernandez-Posada https://orcid.org/0000-0003-3080-1637

C Cochard https://orcid.org/0000-0001-9397-4944

J M Gregg https://orcid.org/0000-0002-6451-7768

R W Whatmore https://orcid.org/0000-0001-9455-5848

M A Carpenter https://orcid.org/0000-0003-2855-0007

## References

[1] Nelmes R J 1974 *J. Phys. C: Solid State Phys.* **7** 3840

[2] Castellanos-Guzman A G 1994 *Proceedings of 1994 IEEE Int. Symp. on Applications of Ferroelectrics* p 6

[3] Schnelle W and Schmid H 2015 *Phys. Rev. B* **91** 184411

[4] Whatmore R W, Herbert J M and Ainger F W 1980 *Phys. Status Solidi a* **61** 73
[5] Schmid H, Genequand P, Pouilly G and Chan P 1980 *Ferroelectrics* **25** 539

[6] Schmid H and Schwarzmüller J 1976 *Ferroelectrics* **10** 283

[7] Sannikov D G 1998 *Ferroelectrics* **219** 177

[8] Schmid H 1994 *Ferroelectrics* **162** 317

[9] McQuaid R G P, Campbell M P, Whatmore R W, Kumar A and Gregg J M 2017 *Nat. Commun.* **8** 15105

[10] Seidel J *et al* 2009 *Nat. Mater.* **8** 229

[11] Farokhipoor S *et al* 2014 *Nature* **515** 379

[12] Van Aert S, Turner S, Delville R, Schryvers D, Van Tendeloo G and Salje E K H 2012 *Adv. Mater.* **24** 523

[13] Rojac T *et al* 2016 *Nat. Mater.* **16** 322

[14] Bednyakov P S, Sturman B I, Sluka T, Tagantsev A K and Yudin P V 2018 *npj Comput. Mater.* **4** 65

[15] Aird A and Salje E K H 1998 *J. Phys.: Condens. Matter* **10** L377

[16] Vasudevan R K *et al* 2012 *Nano Lett.* **12** 5524

[17] Salje E, Devarajan V, Bismayer U and Guimaraes D M C 1983 *J. Phys. C: Solid State Phys.* **16** 5233

[18] Salje E K H 2010 *ChemPhysChem* **11** 940

[19] Catalan G, Seidel J, Ramesh R and Scott J F 2012 *Rev. Mod. Phys.* **84** 119

[20] Sharma P, Zhang Q, Sando D, Lei C H, Liu Y, Li J, Nagarajan V and Seidel J 2017 *Sci. Adv.* **3** e1700512

[21] Dvořák V 1971 *Czech. J. Phys.* **B** **21** 1250

[22] Dvořák V and Petzelt J 1971 *Czech. J. Phys.* **B** **21** 1141

[23] Feng J S, Xu K, Bellaiche L and Xiang H J 2018 *New J. Phys.* **20** 053025

[24] Iliev M N, Hadjiev V G, Íñiguez J and Pascual J 2009 *Acta Phys. Pol. A* **116** 19

[25] Uesu Y, Kobayashi J, Anjoh I and Schmid H 1978 *Ferroelectrics* **20** 167

[26] Kobayashi J, Uesu Y and Kumomi H 1984 *Phase Transit.* **4** 255

[27] Delfino M, Loiacono G M, Smith W A and Gentile P S 1980 *J. Solid State Chem.* **33** 107

[28] Carpenter M A, Salje E K H, Graeme-Barber A, Wruck B, Dove M T and Knight K S 1998 *Am. Mineral.* **83** 2

[29] Meyer G M, Nelmes R J, Thornley F R and Stirling W G 1982 *J. Phys. C: Solid State Phys.* **15** 2851

[30] Lockwood D J 1981 *Ferroelectrics* **36** 443

[31] Moopenn A and Coleman L B 1990 *J. Phys. Chem. Solids* **51** 1099

[32] Carpenter M A, Price G D and Ross N L (ed) 1992 *The Stability of Minerals* (London: Chapman and Hall) p 172

[33] Malcherek T, Kroll H, Schleiter M and Salje E K H 1995 *Phase Transit.* **55** 199

[34] Carpenter M A, Becerro A I and Seifert F 2001 *Am. Mineral.* **86** 348

[35] McKnight R E A, Howard C J and Carpenter M A 2009 *J. Phys.: Condens. Matter* **21** 015901

[36] Whatmore R W, Brierley C J and Ainger F W 1980 *Ferroelectrics* **28** 329

[37] Migliori A and Sarrao J L 1997 *Resonant Ultrasound Spectroscopy: Applications to Physics, Materials Measurements and Nondestructive Evaluation* (New York: Wiley)

[38] Carpenter M A 2015 *J. Phys.: Condens. Matter* **27** 263201

[39] McKnight R E A, Moxon T, Buckley A, Taylor P A, Darling T W and Carpenter M A 2008 *J. Phys.: Condens. Matter* **20** 075229

[40] McKnight R E A, Carpenter M A, Darling T W, Buckley A and Taylor P A 2007 *Am. Mineral.* **92** 1665

[41] Evans D M, Schiemer J A, Schmidt M, Wilhelm H and Carpenter M A 2017 *Phys. Rev.* **B** **95** 094426

[42] Carpenter M A, Bryson J F J, Catalan G, Zhang S J and Donnelly N J 2012 *J. Phys.: Condens. Matter* **24** 045902

[43] Nowick A S and Berry B S 1972 *Anelastic Relaxation in Crystalline Solids* (New York: Academic)

[44] Weller M, Li G Y, Zhang J X, Kê T S and Diehl J 1981 *Acta Metall.* **29** 1047

[45] Carpenter M A, Howard C J, McKnight R E A, Migliori A, Betts J B and Fanelli V R 2010 *Phys. Rev.* **B** **82** 134123

[46] Carpenter M A and Zhang Z 2011 *Geophys. J. Int.* **186** 279

[47] Carpenter M A, Salje E K H and Howard C J 2012 *Phys. Rev.* B **85** 224430

[48] Thomson R I, Jain P, Cheetham A K and Carpenter M A 2012 *Phys. Rev.* B **86** 214304

[49] Zhang Z, Li W, Carpenter M A, Howard C J and Cheetham A K 2014 *CrystEngComm* **17** 370

[50] Xin L, Zhang Z, Carpenter M A, Zhang M, Jin F, Zhang Q, Wang X, Tang W and Lou X 2018 *Adv. Funct. Mater.* **28** 1806013

[51] Carpenter M A, Howard C J, Kennedy B J and Knight K S 2005 *Phys. Rev.* B **72** 024118

[52] Hayward S A *et al* 2005 *Phys. Rev.* B **72** 054110

[53] Carpenter M A, Sinogeikin S V and Bass J D 2010 *J. Phys.: Condens. Matter* **22** 035404

[54] Carpenter M A, Buckley A, Taylor P A and Darling T W 2010 *J. Phys.: Condens. Matter* **22** 035405

[55] Perks N J, Zhang Z, Harrison R J and Carpenter M A 2014 *J. Phys.: Condens. Matter* **26** 505402

[56] Thomson R I, Rawson J M, Howard C J, Turczynski S, Pawlak D A, Lukasiewicz T and Carpenter M A 2010 *Phys. Rev.* B **82** 214111

[57] Spalek L J, Saxena S S, Panagopoulos C, Katsufuji T, Schiemer J A and Carpenter M A 2014 *Phys. Rev.* B **90** 054119

[58] Schiemer J, O'Flynn D, Balakrishnan G and Carpenter M A 2013 *Phys. Rev.* B **88** 054108

[59] Kustov S, Liubimova I and Salje E K H 2018 *Appl. Phys. Lett.* **112** 042902

[60] Catalan G, Saxena A and Planes A (ed) 2014 *Mesoscopic Phenomena in Multifunctional Materials* (Springer Series in Materials Science) vol 198 (Berlin: Springer) p 25

[61] Jia C-L, Mi S-B, Urban K, Vrejoiu I, Alexe M and Hesse D 2008 *Nat. Mater.* **7** 57

[62] Pesquera D, Carpenter M A and Salje E K H 2018 *Phys. Rev. Lett.* **121** 235701

[63] Carpenter M A, Wiltshire E C, Howard C J, Thomson R I, Turczynski S, Pawlak D A and Lukasiewicz T 2010 *Phase Transit.* **83** 703

[64] Carpenter M A, Buckley A, Taylor P A, McKnight R E A and Darling T W 2010 *J. Phys.: Condens. Matter* **22** 035406

[65] Haida M, Kohn K and Kobayashi J 1975 *J. Phys. Soc. Japan* **39** 1625

[66] Rivera J P and Schmid H 1988 *J. Phys. Colloq.* **49** 849

[67] Stokes H T, Hatch D M and Campbell B 2007 *J ISOTROPY Software Suite* (http://stokes.byu.edu/isotropy.html)

[68] Carpenter M A, Salje E K H and Graeme-Barber A 1998 *Eur. J. Mineral.* **10** 621