Magnetic and dielectric order in the kagome-like francisite Cu₃Bi(SeO₃)₂O₂Cl

The first experimental evidence to suggest the presence of a structural distortion in CBSCl was the observation of new phonon modes that appear in the far infrared spectra below $T_{\mathrm{S}}$ = 115 K Miller et al. (2012), a result that was further confirmed by Raman spectra Gnezdilov et al. (2016). This temperature also coincides with a deviation of the magnetic susceptibility from the high temperature Curie-Weiss law, indicating that this structural transition might be correlated with the onset of magnetic correlations Miller et al. (2012).

We first studied the structural distortion using Laue neutron diffraction on CYCLOPS that clearly showed a doubling of the $c$ lattice parameter below $T_{\mathrm{S}}$ (not shown). Using the 2T-1 triple-axis spectrometer at zero energy, we then recorded $\omega$-scans (obtained by rotating the crystal) centered at the scattering vector Q = (2 0 1.5) from 160 K to 5 K. Integration of the neutron counts as a function of temperature are shown in Fig. 2, which indeed reveals the emergence of a Bragg peak below 115 K and is consistent with a structural distortion involving a doubling along the $c$ axis as represented in Fig. 1 (c). Below 25 K, the peak gains further intensity from the magnetic ordering consistent with a ${\bf k}=(0,0,0)$ propagation vector in the new structure with a double $c$ lattice parameter (note that it would correspond to a ${\bf k}=(0,0,1/2)$ propagation vector as for CBSBr without the structural distortion).

To further investigate the nuclear structure below $T_{\mathrm{S}}$, we performed single-crystal neutron diffraction measurements on the D23 diffractometer, collecting Bragg reflections at 150, and 50 K with ${\bf k}=(0,0,0)$ in the orthorhombic cell with a $c$ lattice parameter doubled with respect to the room temperature structure. We carefully inspected the integrated intensities of the Bragg reflections and identified the following reflection conditions in the distorted phase.

Using the SUBGROUPS application of the Bilbao crystallographic server Aroyo et al. (2006a, b, 2011), we identified 10 possible transformations from $Pmmn$ (No. 59) to an orthorhombic space group with a doubling of the $c$ lattice parameter. Of these possible space groups, only centrosymmetic $Pcmn$ (No. 62) and noncentrosymmetric $Pc2_{1}m$ (No. 33) fulfill the above conditions. Notably, the high temperature $Pmmn$ and previously proposed $P2_{1}mn$ space groups fail to meet the observed reflection conditions. We note that possible monoclinic and the triclinic non-centrosymmetric space groups ($Pc$, $Pm$, $P2_{1}$, $P2$ and $P1$) can also be excluded based on the observed reflection conditions. Using the data obtained on D23 at 50 K, we have attempted to fit the integrated intensities in both the $Pcmn$ and $Pc2_{1}m$ space groups using WinPLOTR of the FullProf software package Roisnel and Rodríquez-Carvajal (2001). The quality of the refinements is shown in Fig. 3. The best results are obtained for $Pcmn$ and indicate that displacements of the Cu2 and Cl atoms from $x/a$ = $\frac{1}{4}$ above $T_{\mathrm{S}}$ to $x/a$ = $\frac{1}{4}+\Delta_{\mathrm{Cu2,Cl}}$ below $T_{\mathrm{S}}$, where $\Delta_{\mathrm{Cu2}}$=0.026 and $\Delta_{\mathrm{Cl}}$=0.031 at 50 K, are the defining distortions (see Fig. 1(b)). The refined crystal structure parameters at 2 and 50 K are presented in Table 1 along with the refined parameters at 150 K for the $Pmmn$ structure.

Our results are in agreement with the recent X-ray refinement of the low temperature structure by Prishchenko et al. and discards the polar space group $P2_{1}mn$ proposed initially to account for the additional 16 phonon modes observed below $T_{\mathrm{S}}$ Miller et al. (2012). Actually, a doubling of the $c$ lattice parameter, that was not considered in this previous work, provides a natural explanation to the presence of the new phonon modes as an expansion of the unit cell will increase the number of normal modes it supports. Therefore a loss of centrosymmetry is not required to explain the phonon anomalies.

As shown in the study of Prishchenko et al. Prishchenko et al. (2017), and by the previous literature on francisites including neutron Pregelj et al. (2012) and Raman Gnezdilov et al. (2016) scattering experiments, the structural transition observed in CBSCl is not found in the CBSBr and Cu₃Bi(SeO₃)₂O₂I (CBSI) forms of francisite. Differing only in the types of halogen atoms in the chemical formula, one can reason that it is the size of the halogen atoms in the francisites that is the source of the structural transition in CBSCl. In the francisites, the halogen atom’s position at the center of the kagome hexagons, places it in a potential void where large atomic motion is possible. Unlike Br and I, which are larger atoms than Cl and therefore relatively stable in their positions, Cl will be more free to displace from its position. This, as well as the buckled geometry of the kagome hexagon that places the Cu2 ions out of the kagome plane, promotes an instability of both the Cl and Cu2 ions. Magnetic correlations across the kagome lattice as well as thermal lattice contraction put the Cl and Cu2 ions in energetically unfavorable positions within the $Pmmn$ structure. Therefore, a collective displacement of the Cu2 and Cl ions to a lower energy state occurs. This is consistent with our crystal refinement results, which show antiparallel shifts for Cu2 and Cl along the $a$ crystallographic direction, reducing the Cu2-Cl distance from 3.20 Å  at 150 K to 2.82 Å at 50 K. The new positions are no longer described by the symmetries of the $Pmmn$ space group, rather falling under the $Pcmn$ description with a twice as large $c$ cell parameter. Including this, we also find that as a result of these displacements, small shifts in the O2, O2’ and O3 ions occur. Importantly, this results in a distortion of the SeO₃ triangular pyramid. We note that, while no considerable distortions are observed below our refinement at 50 K, it is the oxygen atoms that have the biggest shifts in position when comparing results at 50 K and 2 K.

To determine the magnetic structure, a further set of neutron Bragg reflections were collected and analyzed at 2 K, below the magnetic transition. While a previous single crystal neutron diffraction refinement of CBSBr has revealed a magnetic structure consistent with the macroscopic magnetic properties of both the Br and Cl compounds Pregelj et al. (2012), no neutron diffraction studies have been reported so far on the Cl compound. One can wonder whether the structural distortion, absent in the Br compound, has an influence on the magnetic order of the Cl compound. In CBSCl the magnetic phase exists within the $Pcmn$ structure which is compatible with a larger set of irreducible representations (IR), with respect to the undistorted case of CBSBr. Through a second order phase transition and for coupled Cu1 and Cu2 spins, one of these IR should correspond to the observed magnetic structure. The 8 possible IRs, determined using BasIreps of the FullProf software package Rodríguez-Carvajal (1993), are reported in Table  2.

We have tested all possibilities and obtained the best refinement for the magnetic structure described by the $\Gamma_{8}$ IR on both Cu1 and Cu2 sites, giving rise to the $Pcm^{\prime}n$ magnetic space group ($Pnm^{\prime}a$ in the standard setting). Contrary to what has been refined for the CBSBr compound, where canting along $a$ of the Cu2 sites is forbidden in the $Pmmn$ structure Pregelj et al. (2012), the $\Gamma_{8}$ IR here allows for a magnetic component along the $a$ axis for both the Cu1 and Cu2 sites. These components are refined by the $u$ and $u^{\prime}$ parameters in Table 2 respectively. In our refinement of CBSCl however, we were unable to experimentally distinguish between a fit where $u$ and $u^{\prime}$ is varied or when they are fixed to $u=u^{\prime}=0$. When $u$ and $u^{\prime}$ are allowed to vary we do indeed refine a structure featuring a very small alternation of the Cu1 spins along $a$ at an angle of 3.8^∘ and $ac$ plane canted Cu2 spins at an angle of 32.6^∘ from $c$ towards $a$. However, such a magnetic structure with finite $u$ and $u^{\prime}$ components disagrees with our mean field calculations and spin wave measurements presented later in the article which suggests no canting along the $a$ axis for either the Cu1 or Cu2 sites. We therefore proceeded using a $\Gamma_{8}$ refined magnetic structure with $u=u^{\prime}=0$ fixed for the Cu1 and Cu2 sites. The quality of the refinement is shown in Fig. 3 (b). The refined magnetic components are then $m_{y}=0.76(3)\mu_{\mathrm{B}}$ and $m_{z}=-0.46(5)\mu_{\mathrm{B}}$ for Cu1 (total magnetic moment of 0.89(4)$\mu_{\mathrm{B}}$) and $m_{z}=m_{tot}=1.03(9)\mu_{\mathrm{B}}$ for Cu2. The magnetic structure is depicted in Fig 4. It is overall similar to the previously reported structure of CBSBr Pregelj et al. (2012): the Cu2 spins are along the $c$ axis and the Cu1 spins are confined to the $bc$ plane alternating along $b$ with a canting of 59^∘ $\pm$4^∘ from $c$ towards $b$. One particular difference is the increased canting of the Cu1 spins with 59^∘ for CBSCl and 50^∘ for CBSBr Pregelj et al. (2012). This difference is possibly a results of the structural distortion in CBSCl, and is discussed later in the context of the magnetic Hamiltonian.

The originality of this magnetic structure, antiferromagnetic stacking of canted spins, was already addressed by Rousochatzakis et al. in relation to the role of magnetic frustration Rousochatzakis et al. (2015). Their first theoretical analysis was aimed at determining the principal magnetic interactions in the francisites using density functional theory. In their study it was found that the nearest neighbor interactions, $J_{1}$ between Cu1 spins along $a$ and $J_{1}^{\prime}$ between Cu1 and Cu2 spins diagonally along $a$-$b$, are dominantly ferromagnetic. A comparably large next nearest-neighbor interaction between the Cu1 spins along $b$, $J_{2}$, is antiferromagnetic and provides the source of frustration. Weak interlayer antiferromagnetic interactions $J_{\perp 1}$ and $J_{\perp 2}$ are also present between the Cu1 and Cu2 sites. The approximate paths of these primary interactions are indicated in Fig 4 (a) and (b). To stabilize the magnetic structure, a DM interaction is introduced between the Cu1 and Cu2 spins with a dominant $a$ axis component that alternates along $b$, depicted by black arrows in Fig 4 (a). Longer range interactions between the Cu1 spins along $a$ and $b$ are determined negligible. The resulting proposed magnetic Hamiltonian is then,

with $\mathbf{D}$ the DM vector and the indexes corresponding to the labeled spins in Fig. 4 (c).

We used a real-space mean-field energy minimization of the spin configuration implemented in the SpinWave software package Petit and Damay (2016) to calculate the magnetic ground state in the $Pcmn$ space group. We note that our simulations produce identical results in the $Pmmn$ framework when adjusting for the different $l$ indices. We therefore began with the above Hamiltonian, using the values $J_{1}$ = -6.55 meV, $J_{1}^{\prime}$ = -5.75 meV, $J_{2}$ = 4.74 meV, $J_{\perp 1}$ = -0.035 meV, $J_{\perp 2}$ = 0.17 meV and $D_{a}$ = 1.04 meV stated in reference Rousochatzakis et al. (2015). Our mean field calculations produce a magnetic structure consistent with the refined $\Gamma_{8}$ structure depicted in Fig. 4. However, using this form of the Hamiltonian we find a smaller Cu1 canting within the $bc$ plane closer to 50^∘ as observed in CBSBr. The Cu2 spins are also not stabilized completely along $c$ with some minor canting along $b$ inconsistent with the representations in Table 2. We note that in reference Rousochatzakis et al. (2015), values for the $D_{b}$ and $D_{c}$ components of the DM interaction have been calculated at $D_{b}\approx D_{c}\approx D_{a}/3$. However, in our simulations we found that while the $D_{b}$ and $D_{c}$ components were required to stabilize the magnetic structure, their values could be much smaller without affecting the minimization. In fact at the values stated in reference Rousochatzakis et al. (2015), the $D_{b}$ and $D_{c}$ components had the effect of distorting the calculated spin wave dispersions. We found that a more appropriate relation is given by $D_{b}\approx D_{c}\approx D_{a}/10$. Finally, starting from the ground state and minimizing the energy, we are able to simulate the spin wave spectrum shown in Fig. 5 between 0–6.5 meV for $\mathbf{Q}=(0,k,1)$ (a) and $\mathbf{Q}=(0,1,l)$ (b) in good agreement with the simulations of reference Rousochatzakis et al. (2015).

We have compared this model to our inelastic neutron scattering measurements. Two series of spin wave spectra were measured in the ($b^{*}$, $c^{*}$) scattering plane at 2 K, from 0–8 meV with energy steps of 0.2 meV, along $\mathbf{Q}=(0,k,2)$ and $\mathbf{Q}=(0,1,l)$ with $\Delta Q_{h}=0.125$ and $\Delta Q_{l}=0.5$ intervals. Selected spectra are shown in Fig. 6 (a-d). They are combined into the dispersion maps displayed in Figs. 5 (c) and (d). The spin wave branches below 8 meV are dispersive along $k$ and gapped. At the zone centers, the resulting excitation spectra feature peaks at 1.57 meV and 4.34 meV (see Fig. 6(a)), consistent with the results of THz Miller et al. (2012) and Raman Gnezdilov et al. (2016) spectroscopy. This disagrees however with the spin wave analysis of Rousochatzakis et al. Rousochatzakis et al. (2015), who found gapless excitations at the zone centers, indicating some inaccuracies in the Hamiltonian of equation 1. The spin waves are dispersionless along $c$ within the energy resolution given by the experimental conditions (see Figs. 6(a) and  5(d)). This highlights the 2-dimensional nature of the kagome planes due to the weakness of the interlayer coupling. When heating to 15 K, both low energy modes at the zone centers shift toward lower frequencies, consistent with temperature dependent Raman scattering spectra Gnezdilov et al. (2016). Both modes are completely absent at 30 K and higher temperatures, hence above $T_{\mathrm{N}}$ = 25 K, confirming their magnetic origin.

Many of the key features of these measured spin waves, including the dispersion and the spectral weight of the excitations, are captured in the simulations produced by the Hamiltonian of Rousochatzakis et al. Rousochatzakis et al. (2015) (equation 1). However, as mentioned previously, the spin gap is absent. To reproduce this gap we have attempted many modifications to the Hamiltonian including adjusting the interaction strengths, altering the DM interaction vector and introducing anisotropy. We also considered the consequences on the magnetic Hamiltonian of the structural distortion occurring at 115 K from the $Pmmn$ to $Pcmn$ space groups with a doubling of the unit cell along $c$. From our neutron scattering refinements, we find that the Cu1-Cu2 bond distance shifts from a symmetric 3.25 Å at 150 K to either 3.19 Å or 3.32 Å at 2 K. This leads to a splitting of the $J_{1}^{\prime}$ interactions into two interactions, $J_{11}^{\prime}$ for the shorter bonds and $J_{12}^{\prime}$ for the longer bonds. This is also expected to split the DM interactions. When assessing the refined Hamiltonian, comparisons were made between the mean-field stabilized magnetic structure and our refined structure of Fig 4. The simulated spin waves were then compared with the measurements, checking peak positions and intensity over the energy-$Q$ space probed in our experiments. Finally, consideration was also given to whether the adjustments were physically acceptable or not. However, we found little change from these modifications to the simulated magnetic properties, including no induced gap in the spin wave excitations and limited increased canting of Cu1 spins.

More successful attempts to reproduce the observed magnetic structure and the associated spin wave dispersion were actually achieved by introducing an adapted anisotropic exchange. The symmetric portion of the anisotropic exchange, that adds to the isotropic exchange interactions of equation 1, is implemented into the Hamiltonian with the following term,

where $\mathbf{\Gamma}_{i,j}$ is a second-rank symmetric tensor. The tensors $\mathbf{\Gamma}_{1,2}$, $\mathbf{\Gamma}_{5,2}$ and $\mathbf{\Gamma}_{3,2}$ then represent the symmetric anisotropic exchange for the bonds between the sites labeled in Fig.4 (c) (as well as all other symmetrically equivalent bonds). In their calculations, Rousochatzakis et al. found that the contribution of this symmetric anisotropic exchange was negligible Rousochatzakis et al. (2015). However, in our analysis we have found that the combination of small positive (antiferromagnetic) off-diagonal $\Gamma^{yz}_{1,2}$ and $\Gamma^{zy}_{1,2}$ components and larger negative components for $\Gamma^{xx}_{4,2}$ and $\Gamma^{zz}_{4,2}$ is able to open the gap in the spin wave spectrum. Further adjustment of the gap energy and of the global dispersion, as well as the canting of the Cu1 spins, is achieved by including a small component for $\Gamma^{xx}_{5,2}$ and by slightly changing the other terms of the Hamiltonian. Finally, we find an adequate solution by only a slight increase in the DM interaction accompanied by an adjustment of $J_{2}/\left|J_{1}^{\prime}\right|\approx 1$, yielding the following modified magnetic Hamiltonian with $J_{1}$ = -6.55 meV, $J_{1}^{\prime}$ = -5.75 meV, $J_{2}$ = 5.7 meV, $J_{\perp 1}$ = -0.035 meV, $J_{\perp 2}$ = 0.17 meV and $D_{a}$ = 1.3 meV, with $D_{b}=D_{c}=D_{a}/10$, and the following symmetric portion of the anisotropic exchange with values represented in meV:

The terms that are significant for generating the spin wave gap are underlined. Note that the anisotropic exchange tensors are essentially equivalent to those calculated in Ref. Rousochatzakis et al. (2015) with the exception of the relatively large components we have found for $\mathbf{\Gamma}_{4,2}$. We note that there is a high degree of uncertainty associated to the refined values in $\mathbf{\Gamma}_{4,2}$. Taking the presence of a spin wave gap and a $\sim 3$ meV splitting of the two spin wave branches at $k=1$ as refinement limits, we find an approximate upper boundary for the uncertainties of these components to be $\mathbf{\Gamma}^{xx}_{4,2}=-0.7\pm 0.4$ and $\mathbf{\Gamma}^{zz}_{4,2}=-1.7\pm 1.0$.

By incorporating this anisotropic exchange into the Hamiltonian, we are thus able to reproduce the magnetic structure with no instability along $b$ for the Cu2 spins. The calculated canting of the Cu1 spins is 59^∘ from $c$ to $a$ in close agreement with our single crystal neutron diffraction refinement. No canting is found along the $a$ direction. The simulated spin wave dispersions using the modified Hamiltonian for $\mathbf{Q}=(0,k,2)$ and $\mathbf{Q}=(0,1,l)$ are shown in Figs. 5 (e) and (f) respectively. Here we find excellent agreement between experiment (Figs. 5 (c), (d)) and simulation (Figs. 5 (e), (f)) when comparing the intensity and dispersions along $k$ and $l$. Our INS experiments thus reveal that a next-nearest-neighbor exchange between the Cu1 sites in the $ac$ plane includes a non-negligible anisotropic exchange part. This then induces the gapped spin wave spectrum we observe.

Further agreement between experiment and simulation is obtained for the ($a^{*}$, $c^{*}$) scattering plane. As is predicted in our simulations, the scattering intensity of the spin waves for the ($a^{*}$, $c^{*}$) scattering plane are reduced by a factor of $\sim$3 when compared to the ($b^{*}$, $c^{*}$) scattering plane, rendering it more difficult to map the dispersion in this orientation. Nevertheless, by combining measurements performed on several spectrometers, we were able to extract peak positions across a broad span of $h$ between 0 and 20 meV. The resulting experimental and simulated dispersion is shown in Fig. 7 (a). The extracted peak positions, plotted on top of the simulation using the modified Hamiltonian, show a good agreement at energies below 8 meV. Fig. 7 (b) reveals the predicted dispersion of higher energy spin waves in the ($b^{*}$, $c^{*}$) scattering plane including extracted peak positions above 8 meV. The energies of transverse optical phonons observed by infrared spectroscopy Miller et al. (2012) are indicated by dashed gray lines plotted over the simulations.

Concerning the higher energy spin wave branches in both the ($a^{*}$, $c^{*}$) and ($b^{*}$, $c^{*}$) scattering planes, no significant magnetic contribution could be detected in the INS experiments when comparing spectra above and below $T_{\mathrm{N}}$. In fact, while many of the higher energy peaks (above 8 meV) measured in our INS experiments are aligned with predicted spin wave branches, they are also in close proximity to the previously observed infrared optical phonons. The presence of spin-wave–phonon hybridization has been previously established in CBSCl Gnezdilov et al. (2016). It is therefore possible that the higher energy spin wave branches are either masked or hybridized with phonon signatures observed in Fig.7. Further experiments incorporating polarized neutron analysis will be required to investigate this possibility.

To further study the phase diagram of CBSCl, we have looked at possible signatures of the structural and magnetic transitions in the dielectric properties. Electrical capacitative measurements were used to determine the relative permittivity of CBSCl shown in Figs. 8 (a) and (b) for a single crystal along the $c$ axis and for a powder pellet sample respectively. In Fig. 8 (a), a dielectric anomaly is observed at $\sim$120 K, coinciding with the structural transition. The anomaly is continuous, which is characteristic of a 2^nd-order phase transition. At the dielectric anomaly, we observe no frequency dependence, thus indicating long-range polar correlations, nor thermal hysteresis suggesting no kinetic effects. Since we have shown that the structural transition involves a distortion from the centrosymmetric $Pmmn$ space group to the centrosymetric $Pcmn$ space group, we rule out the formation of long-range polar order into a ferroelectric (FE) phase and turn to the possibility of an antiferroelectric (AFE) low temperature phase.

Antiferroelectrics are a somewhat loosely defined class of non-polar materials. They are typically distinguished from standard non-polar systems at the microscopic scale by considering the presence of non-collinear dipoles mutually canceling Tolédano and Guennou (2016). They tend to form following a distortion from one symmetric phase to a less symmetric phase accompanied by an anomaly in the dielectric permittivity related to the formation and/or organization of local dipole moments in an AF arrangement. This can also involve a doubling of the cell parameters when the associated soft phonon modes feature two groups of displacements, equal and opposite in magnitude Rabe (2013); Devonshire (1954). The basic theoretical description of these AFEs is therefore one with two (or multiples of two) equivalent and opposing polarized sublattices Kittel (1951). The application of a sufficiently large electric field along the dipolar axis is expected to invert the polarization of the opposing sublattice, forming an induced ferroelectric state, and leading to a characteristic hysteresis loop in the electrical polarization as a function of the electric field Tolédano and Guennou (2016); Rabe (2013). The $Pcmn$ phase in CBSCl is consistent with the transition at $T_{\mathrm{S}}$ from a paralelectric (PE) phase to an AFE state as described above. The primary distortion involves displacements of the Cu2 and Cl ions in opposing directions along $a$ giving rise to a doubling of the $c$ cell parameter. Note that the isotropic displacement parameter of the Cl^- ion obtained from the neutron refinement (see Table 1) is relatively large, which can indicate a high degree of mobility. This distortion is accompanied by the anomaly in the dielectric permittivity and provides clear candidates for opposing dipoles lying primarily within the $ab$ plane.

The thin planar geometry of the single crystal samples confines electrical measurements along the $c$ axis. To probe the dielectric response in the $ab$ plane, we look at measurements on the powder pellet sample, averaged over all crystallographic directions. As evident in Fig. 8 (b), the same dielectric anomaly at 120 K is observed in the powder pellet sample signifying the formation of the AFE phase. Again, no frequency dependence is present. However, in contrast to the single crystal results, a second anomaly is observed at higher temperatures around 230 K displaying thermal hysteresis as well as a small frequency dependence. The origin of this feature could arise from different dynamics for highly mobile electric charges present in the material, for instance the Cl^- ions, which are impeded by the many grain boundaries of the pellet sample but not for the single crystal.

To provide further information regarding the high temperature dispersive anomaly of Fig. 8 (b) and the nature of the PE-AFE transition in general, we performed additional electric polarization measurements. The temperature dependent depolarizing current measurements used to calculate the electrical polarization are shown in Fig. 9 for the single crystal sample. A large current is detected above $T_{\mathrm{S}}$ forming a peak around 160 K. It follows a linear dependence on the poling bias and is flipped when the polarization of the bias is reversed. It is related to the dispersive dielectric anomaly discussed above, which appears at higher temperatures but shifts to lower temperatures at lower frequencies, consistent with the occurrence at 160 K of the current peak, which can be considered a DC measurement. The released current is detected in both the single crystal and powder pellet sample (not shown), which suggests that the process is intrinsic to the material and does not depend on the sample form. This depolarizing current may be the result of trapped charges that are freed when microscopic AFE domains relax in the unstable thermal region above $T_{\mathrm{S}}$. Their contributions to the pyroelectric current are removed when considering calculations of the electrical polarization.

The pyroelectric current, measured for an electric field applied along the $c$ axis on a single-crystal, is shown in the inset of figure 9 for the low temperature phases. Two small features are detected, one at the structural transition ($T_{\mathrm{S}}$) and one at the magnetic transition ($T_{\mathrm{N}}$). The resulting electrical polarization over this region is shown in Fig. 10. A small electrical polarization is detected below $T_{\mathrm{S}}$, increasing at lower temperature, with an additional small contribution occurring below $T_{\mathrm{N}}$. At electric fields up to $\pm$200 kV/m, neither the contribution below $T_{\mathrm{S}}$ nor the contribution below $T_{\mathrm{N}}$ are affected by the poling bias. At higher electric fields in the range of $\pm$1000 kV/m, the polarization is manipulated below $T_{\mathrm{S}}$, switching polarity between the application of positive and negative poling fields. However, the contribution to the polarization below $T_{\mathrm{N}}$ remains unaffected, even at the highest poling fields achievable in our experimental setup. The polarization obtained for the powder pellet sample is shown in Fig. 10 (b). The polarization increases for decreasing temperature, and is significantly larger than that detected for the single crystal along the $c$ axis under equivalent poling fields. Under small poling fields a negligible polarization is observed, likely only displaying a contribution from the residual polarization detected along the $c$ axis. Under larger poling fields, the polarization is increased and is switchable with the polarity of the field. Notably, the polarization contribution below $T_{\mathrm{N}}$ detected along the $c$ axis is not observed in the powder pellet measurement.

In order to interpret these measurements, we recall that the structural distortion at 115 K involves a transformation between two centrosymmetric space groups, $Pmmn$ and $Pcmn$, that do not support macroscopic polar moments. Actually, the polarization switching in the polycrystalline sample and at high poling fields in the single-crystal is consistent with an induced polarization in an AFE phase, and with the notion of antiparallel dipoles in the $ab$ plane due to the displacements of the Cu2 and Cl ions along $a$.

Another perplexing result is the presence of the weak non-switchable polarization at $T_{\mathrm{N}}$ with no signature in the dielectric permittivity. It is characteristic of some form of magnetoelectric coupling, and only detected in the single-crystal. The fact that is is not detected in the powder sample suggests that the magnetoelectric polarization is indeed along the $c$ axis. The magnetic point group $mm^{\prime}m$ that describes the magnetic phase of CBSCl is compatible with a linear magnetoelectric coupling, $\alpha^{ME}_{ij}E_{i}H_{j}$, with $E_{i}$ and $H_{j}$ the components of the electric and magnetic fields respectively and $\alpha^{ME}_{ij}$ the components of the magnetoelectric tensor of the form Borovik-Romanov and Grimmer (2006):

The off-diagonal terms indicate that the coupling is maximized when the electric and magnetic fields are perpendicular. Its form also allows ferrotoroidicity as the toroidic moment is proportional to the antisymmetric part of the magnetoelectric tensor Gorbatsevich et al. (1983). However in the absence of a uniform magnetic field, it should not lead to a net macroscopic polarization. This anomalous electric polarization we observe may rather be associated to symmetry breaking of local magnetic fields generated at defect or domain wall boundaries. A recent study of the multiferroic behavior in CBSCl has indeed shown that an anisotropic spin-flip transition induced by magnetic fields above $\sim$0.8 T leads to a breaking of the magnetic two-fold symmetry in the $bc$ plane Wu et al. (2017). The resulting phase demonstrates ferroelectric behavior below $T_{\mathrm{N}}$.