Raman imaging of ferroelastically configurable Jahn–Teller domains in LaMnO₃

Abstract

The Jahn–Teller (JT) effect, through geometric deformation of molecules or local ionic lattices, lowers the overall energy of the system by removing electron degeneracy from partially occupied orbitals. Crystal symmetry lowered by JT distortion inevitably creates multiple variants of elastic and orbital-anisotropic states. Visualization and control of the domain/wall textures create a cornerstone to understand various correlated phenomena and explore wall properties. Here, we report the real-space observation of JT phonon and orbiton-related domains in a LaMnO₃ thin film using confocal Raman spectromicroscopy. The characteristic symmetries of the JT-originated Raman modes allow us to detect and visualize the local population and orientation of the JT planes. Combined with a crystal structural analysis, we find that the formation of ferroelastic domains with W or W’ walls provides the basic framework for understanding JT domain textures. Furthermore, we demonstrate the JT domains can be manipulated by applying local external stress. Our findings provide a useful pathway for mechanically-tunable orbitronic applications.

Introduction

Transition metal oxides exhibit a plethora of functional properties and quantum phenomena, including high-temperature superconductivity, colossal magnetoresistance, and electrochemical catalytic activity via strong inter-couplings among spin, charge, and lattice degrees of freedom^1,2,3,4,5,6. Double or superexchange interactions were proposed to account for the interplay of magnetism and transport properties^7,8, and the basic principles strongly rely on the presence of charge ordering^9,10, Jahn–Teller (JT) distortion^11,12,13, and orbital ordering (OO)^14,15,16,17. The JT distortion involves the partial removal of the d-electron degeneracy of a transition metal cation in a ligand field, resulting in a reduction in lattice symmetry. Since a system with orbital degeneracy is electronically unstable, the surrounding anions are deformed to remove the degeneracy. Cooperative manifestation of such local anion distortions leads to a long-range coherence, known as the cooperative JT distortion, which is necessarily accompanied by a long-range order of anisotropic orbitals, that is, OO.

The cooperative JT distortion is largely responsible for the complex relationship between the electronic and magnetic properties of the OO states^17,18,19. For example, epitaxial deformation of La_xCa_1-xMnO₃ thin films promotes JT distortion, suppressing ferromagnetic properties and stabilizing the antiferromagnetic insulating phase^20,21. A ferromagnetic phenomenon due to a self-doping mechanism was also observed at the polar-nonpolar interface of LaMnO₃(LMO)/SrTiO₃(STO) heterostructures^22,23. Moreover, doped manganites undergo an inherent phase competition, causing the quantum ground states to fluctuate spatially. Accordingly, the appearance of domains, walls, and spatially non-uniform textures is an intriguing feature that strongly influences electrical and magnetic properties^9,24,25. Similarly, the metal-insulator transition observed in LMO under high pressure is also understood based on a phase-separated configuration consisting of spatial regions of JT distorted and undistorted MnO₆ octahedra^26,27. Theoretical analysis predicts a conduction enhancement at the anti-phase JT domain wall in an undoped manganite²⁸.

To unravel the role of JT distortion and OO in strongly correlated functionalities and implement various proposed phenomena^11,25,28, the observations of cooperative JT distortions and concomitant OO states are necessary. The fundamental properties have been successfully explored using resonant x-ray scattering^14,15,16 or through detection of characteristic phonon modes using Raman spectroscopy^26,29. However, real-space access to the domain and wall textures of JT distortion along with nanoscale engineering and domain texture control is still difficult to achieve.

We employ confocal Raman spectroscopy, which is one of the far-reaching imaging systems capable of collecting high-quality multispectral images with a submicron lateral resolution and is widely used in a variety of research fields^30,31,32. Its applicability is attributed to technological advances in high-efficiency laser sources, low-noise detectors, effective Rayleigh filters, and high-throughput optics (Fig. 1a)³³. As a fast, convenient, and non-destructive technique, Raman spectroscopy has been used to directly observe the vibrational/rotational states of given crystal symmetries^26,29,34 and to confirm the elementary excitations of an order parameter such as magnon or orbiton^35,36,37.

Fig. 1: Polarized Raman spectra and maps of a 150-nm-thick LMO thin film.

a Schematic of the experimental setup for the confocal Raman spectroscope with backscattering geometry. b Schematic diagrams of a JT distortion and the resultant energy shift of Mn e_g levels, and the concomitant OO of singly occupied e_g orbitals with the cooperative octahedral deformation within the ab plane. The ab JT planes are stacked along c axis and spaced from each other with a relatively smaller pseudocubic lattice parameter. c Raman spectrum measured in (yy) polarization geometry. In the inset, the anion vibrational patterns of two octahedra are depicted for A_g and B_2g modes. d $C_2^x$ symmetry operation (two-fold rotation around a-axis) on the vibrating LMO cells under the B_2g symmetry results in the opposite phase of the vibration direction, denoted by the character (χ) of −1. e AFM image. f Raman map of the A_g mode at ~490 cm⁻¹. g Raman map of the B_2g mode at ~610 cm⁻¹. The areas marked in e-g with solid and dashed circles indicate the identical regions, respectively. The dark area in the lower right corner of the images is an intentionally burned area using a strong laser intensity. Scale bar indicates 2 μm.

The perovskite LMO (space group Pnma) is a canonical system for OO harmonized with a strong JT distortion, by the hybridization of the Mn 3d and O 2p orbitals, wherein four electrons occupy Mn³⁺ d-orbitals with a high spin configuration (t_2g³e_g¹)³⁸. Two-fold degenerate e_g orbitals undergo the JT distortion (T_JT ~750 K)³⁹, which lifts the degeneracy and lowers the electronic energy by single electron occupation (Fig. 1b). The OO state is established by cooperatively distorted oxygen ion cages connected to each other via corner sharing, resulting in the alternative occupation of the $d_{3x^2 - r^2}$/$d_{3y^2 - r^2}$orbitals in the orthorhombic ab plane^14,15 as schematically shown in Fig. 1b. A-type antiferromagnetic order is established via the super-exchange interaction between ordered inter-site e_g electrons at T_N ~140 K³⁸. In the bulk LMO with a D_2h point group symmetry, group theoretical analysis leads to 24 Raman active modes (7 A_g + 5 B_1g + 7 B_2g + 5 B_3g) among 60 Γ-point phonon modes²⁹. According to lattice dynamical calculations²⁹, it was shown that the phonon B_2g (A_g) mode near ~610 cm^‒1 (~490 cm^‒1) corresponds to the symmetric (antisymmetric) stretching of basal oxygens of the MnO₆ as shown in Fig. 1c. The phonon mode with A_g symmetry is associated with both the JT distortion and the octahedral rotation, but the phonon mode with B_2g symmetry originates purely from the cooperative JT distortion (Fig. 1c). The crucial effect of the JT distortion on the symmetric stretching mode can be identified from the symmetry analysis based on the character table of D_2h point group (see Fig. 1d and Supplementary Table. 1). The symmetry operation of $C_2^x$ on a lattice with antiferrodistortive JT distortion results in a reversal of the stretching directions confirming the character (χ) of −1, but it becomes +1 in the case without the JT distortion.

In this paper, we implement real-space imaging of the cooperative-JT-distortion and OO domains in an epitaxial LMO thin film using polarized Raman spectroscopy. The interrelationships between dense striped JT domains, surface topographic features, and orthorhombic crystal structures are carefully investigated in an integrated manner. We discuss the JT domain construction in terms of ferroelasticity and then demonstrate the elastic manipulation of the JT domains by applying local external stresses. This discovery provides a useful route for orbital engineering of complex oxides.

Results and discussion

Domain structure in a LaMnO₃ thin film

We prepared a ~150-nm-thick LMO thin film on a (001)_C-oriented Nb-doped (0.5 wt%) STO substrate (where the subscript ‘_C’ stands for the crystal axes of the STO substrate in the cubic index) by pulsed laser deposition under a high-temperature growth condition that maintains OO even in thin films¹⁵. Figure 1c shows the polarized Raman spectrum that is measured in the (yy) polarization geometry (i.e., both incident and scattered lights have a linear polarization along y axis) and averaged over the whole region of a topographic image in Fig. 1e. Here, the x, y, and z directions are taken to be the directions of the crystallographic orientations of the substrate, i.e., [100]_C, [010]_C, and [001]_C, respectively. The Raman spectrum shows peak features that have been observed in a previous report²⁹, and the symmetries of the peaks can be assigned as aforementioned. Raman spectral peaks in the higher energy regime (>800 cm⁻¹) have been suggested to be related to the elementary excitations of OO³⁵, but the possible involvement of a lattice mode has not been completely ruled out^{40,41,42,43,44}.

Figure 1e–g display an atomic force microscope (AFM) image and Raman maps of A_g (~490 cm⁻¹) and B_2g (~610 cm⁻¹) modes over an identical region. As clearly seen in Fig. 1e, the domain structures of the LMO film are composed of topographically flat regions and modulating regions along the {110}_C directions. Considering the lattice parameter relationship between the bulk LMO (a_bulk = 3.990 Å, b_bulk = 3.990 Å, c_bulk = 3.846 Å in terms of a pseudocubic unit cell)³⁹ and the STO (a_STO = b_STO = c_STO = 3.905 Å) substrate⁴⁵, i.e., a_bulk > a_STO > c_bulk, we presume that the ab JT planes (Fig. 1b) are most likely vertical to the substrate and a mosaic mixture of the vertical planes with four-variant orientations is responsible for the topographic feature through reducing misfit strain energy with the substrate.

The Raman map of the A_g mode at ~490 cm⁻¹ (Fig. 1f) has a nearly uniform distribution. In contrast, the Raman map of the B_2g mode at ~610 cm⁻¹ (Fig. 1g) exhibits clear stripe domains with alternating higher and lower intensities. It is worth noting that the flat region (marked in the solid circle) has clear stripe domains in the Raman map of the B_2g mode, while the topographically modulated region (marked in the dashed circle) does not have distinctive stripes. We will analyze the Raman signals to clarify the correlation with the topographic feature in the following sections.

Angle-resolved Raman spectroscopy

Because of the possible multiple orientations of the JT planes (ab planes) in epitaxial thin films¹⁵, we first assess the symmetry of Raman intensity quantitatively concerning the possible geometries of the JT plane based on the group theory. Based on the symmetry analysis, the Raman tensor R of the B_2g and A_g modes can be expressed in the light polarization space as

$$R_{\mathrm{O}}(A_g) = \left( {\begin{array}{*{20}{c}} {\alpha _{xx}} & 0 & 0 \\ 0 & {\alpha _{yy}} & 0 \\ 0 & 0 & {\alpha _{zz}} \end{array}} \right),$$

(1)

and

$$R_{\mathrm{O}}(B_{2g}) = \left( {\begin{array}{*{20}{c}} 0 & {\alpha _{xy}} & 0 \\ {\alpha _{yx}} & 0 & 0 \\ 0 & 0 & 0 \end{array}} \right).$$

(2)

Here, each tensor element α_ij (i, j = x, y, z) interrelates the incident light intensity of electric polarization (j) with the Raman scattered radiation of polarization (i). When we set the bases of the matrix as the principal axes of the orthorhombic LMO unit cell (note the subscript ‘_O’ represents an orthorhombic notation as described in Fig. 1b), the two Raman matrices can be purely diagonal or off-diagonal, and we approximate that the diagonal components of R_O(A_g) have an identical value, i.e., α_xx ≈ α_yy ≈ α_zz ≡ α considering the perovskite unit cell has a shape of pseudocube. Accordingly, the B_2g mode can play a decisive role in determining the possible orientations of the JT plane. Using Euler angles composed by elemental rotations to transform the lab frame (L) described in Fig. 1a into the orthorhombic frame (O), $S_{{\mathrm{L}} \to {\mathrm{O}}}$, we can find the Raman tensor in the lab frame for the given orientation of orbital ordering planes,

$$R_{\mathrm{L}} = S_{{\mathrm{L}} \to {\mathrm{O}}}^\dagger R_{\mathrm{O}}S_{{\mathrm{L}} \to {\mathrm{O}}}$$

(3)

and Raman intensity can be obtained using the following equation

$$I \propto \left(\mathop {\sum}\limits_{i,j} {e_i\alpha _{ij}e_j} \right)^2.$$

(4)

Here, e_i and e_j represent the polarization vector components of the scattered and incident light, respectively. Based on our previous report¹⁵, we consider the symmetry of the Raman modes for the three representative orientations of the JT plane as candidates. Figure 2 shows the JT planes parallel to the (100)_C, (010)_C, and (001)_C planes (Fig. 2a–c) and the resultant angle-resolved Raman intensities for two different polarization geometries (Fig. 2d–g). Here, the azimuthal angle (φ) in the polar plots spans from the substrate crystallographic [100]_C direction to the polarization direction of the normal-incident light (Fig. 2). Figure 2d, e display the Raman intensities of the B_2g mode for the given orientations of the JT planes. If the JT planes are perpendicular to the substrate surface, the Raman intensities in the parallel polarization geometry become

$$I_\parallel (B_{2g}) \propto \left\{ {\begin{array}{*{20}{c}} {{\mathrm{sin}}^4\varphi \;\;{\mathrm{for }}\,(100)_{\mathrm{C}}\;{\mathrm{JT}}\,{\mathrm{plane}}} \\ {{\mathrm{cos}}^4\varphi \,{\mathrm{for }}\,(010)_{{\mathrm{C}}\,}\;{\mathrm{JT}}\,{\mathrm{plane}}} \end{array}{\mathrm{.}}} \right.$$

(5)

Fig. 2: Theoretical analysis of the Raman intensities.

Schematic diagrams describing the possible orientations of JT planes parallel to substrate (100)_C (a), (010)_C (b), and (001)_C (c). Polar plots of the B_2g modes in parallel (d) and cross (e) polarization configurations. Polar plots of the A_g modes in parallel (f) and cross (g) polarization configurations.

With two-fold symmetry, the maximum intensity is obtained when the light polarization is aligned along the JT planes. Meanwhile, the four-fold symmetry can be obtained when the JT planes are parallel to the substrate surface as

$$I_\parallel (B_{2g}) \propto \left| {\cos 2\varphi } \right|^2\,{\mathrm{for}}\,(001)_{\mathrm{C}}\,{\mathrm{JT}}\;{\mathrm{plane}}.$$

(6)

In the cross polarization geometry (the polarization of incident light is perpendicular to that of scattered light), all of the cases have four-fold symmetry, but the intensity is significantly reduced when the JT plane is perpendicular to the substrate, i.e.,

$$I_ \bot (B_{2g}) \propto \left\{ \begin{array}{l}{\textstyle{1 \over 4}}\left| {\sin 2\varphi } \right|^2\,{\mathrm{for}}\,(100)_{\mathrm{C}}\,{\mathrm{or}}\,(010)_{\mathrm{C}}\,{\mathrm{JT plane}}\\ \left| {\sin 2\varphi } \right|^2\,{\mathrm{for}}\,(001)_{\mathrm{C}}\,{\mathrm{ JT plane}}\end{array} \right..$$

(7)

Meanwhile, the Raman intensity $I_\parallel (A_g)$ or $I_ \bot (A_g)$ of the A_g mode shows no dependence of the azimuthal angle. $I_\parallel (A_g)$ has a non-zero value in the parallel polarization geometry (Fig. 2f), however, $I_ \bot (A_g)$ vanishes in the cross-polarization geometry (Fig. 2g) regardless of the orientations of the JT planes. This is the reason why the Raman map in the B_2g mode has a clear color contrast, in contrast to the uniform A_g mode map in Fig. 1f, g.

Figure 3 and Supplementary Fig. 1 show the experimental results from the angle-resolved Raman spectroscopy of the B_2g mode for three representative positions, which are indicated by the Raman maps represented above the polar plots. For the quantitative comparison, the intensities of the B_2g mode were normalized by those of the A_g mode. Figure 3a–c (Fig. 3d–f) show the polar plots obtained in the parallel (cross) polarization geometry. The scattered data from 0° to 180° are experimental data (solid circles) and the data from 180° to 360° (empty circles) are a repetition of that from 0° to 180°, to provide a better comparison between the polar diagram and our fitting curves.

Fig. 3: Experimental results of the Raman intensities.

a–f Polar plots of the Raman intensities of the B_2g mode normalized by that of the A_g mode at three different positions, marked in the Raman maps above the corresponding graphs. The Raman maps were obtained at ~ 610 cm⁻¹ in the (xx) polarization configuration. The polar plots measured in the parallel (a–c) and cross (d–f) polarization geometries are displayed. The solid circles from 0° to 180° are experimental data, and the empty circles from 180° to 360° are the repetition of the data from 0° to 180° in order to construct the polar diagram. g The fitting functions consist of two planes perpendicular to each other. The blue (red)-boxed terms indicate (010)_C ((100)_C)-oriented JT planes, as displayed in the schematics above the formula. h The definition of the director, which is parallel to the JT planes. Scale bar in the Raman map is 1 μm.

In the bright region of the Raman map, the intensity has a two-fold symmetry, with a spindle shape where the intensity is maximized when the light polarization is aligned along the [100]_C direction (Fig. 3a). In the dark region, the polar plot in Fig. 3b has a similar shape to the one in Fig. 3a but the graph is rotated by 90°. The two-fold symmetries suggest that the majority of the signals are from the (100)_C or (010)_C JT planes, respectively. Meanwhile, the area with the intermediate intensity has four-fold symmetry (Fig. 3c). The four-fold symmetry seems to come from the (001)_C JT planes or from the composite structures of the (100)_C and (010)_C JT planes. The former case is excluded by the results of Fig. 3d–f, since the polar plots for all positions in the cross polarization case have the same symmetry and similar Raman intensities. Therefore, the JT planes of our LMO film prefer to be perpendicular to the substrate, and the symmetry of the polar plots can be attributed to the relative domain population of the (100)_C and (010)_C JT planes. Therefore, we set our fitting model for the polar plots based on the Eq. (5) with possible in-plane rotation $\bar \varphi$ of the JT vertical planes.

$$f = \left\{ \begin{array}{l}f_{||} = a_{||}{\mathrm{ cos}}^4(\varphi - \bar \varphi ) + b_{||}{\mathrm{ sin}}^4(\varphi - \bar \varphi )\quad \quad {\mathrm{for}}\, {\mathrm{paralell}}\, {\mathrm{polarization}}\, {\mathrm{geometry}}\\ f_ \bot = a_ \bot {\mathrm{ sin}}^22(\varphi - \bar \varphi ) + b_ \bot {\mathrm{ sin}}^22(\varphi - \bar \varphi )\quad {\mathrm{for}}\, {\mathrm{cross}}\, {\mathrm{polarization}}\, {\mathrm{geometry }}{\mathrm{.}}\end{array} \right.$$

(8)

The first (second) term in Eq. (8) describes the azimuthal-polarization-angle dependence for the (010)_C ((100)_C) oriented JT planes, as illustrated in Fig. 3g. Our model fits all of the experimental data well. Since the Raman intensity is proportional to the domain population³⁴, the directional preference of the JT plane within a specified position can be determined. In this regard, we define a local director, which is a vector-like quantity with two equivalent opposite directions (headless vector), to understand the spatial distribution described in Fig. 3h. Such a director is widely used as an order parameter in nematic liquid crystals, where the local preferential order of the molecules affects the phase of the liquid crystals. In our system, the director describes the preferred direction of the JT planes, and the length of the director indicates the relative population of two perpendicular domains.

Real-space director field

To construct the real-space textures of the directors associated with the local orientation of the JT planes, we first established Raman maps over the same region for various polarization angles in the parallel polarization geometry. To align the multiple Raman images pixel by pixel, we used a modified image registration algorithm for optimizing the conversion matrices among images⁴⁶. Based on our analysis method, introduced in the aforementioned section, the Raman maps of the B_2g mode were normalized by those of the A_g mode, as shown in Supplementary Fig. 2. We carried out fitting using Eq. (8) to find the coefficients per each position and mapped out the directors spatially.

We successfully visualized the distribution of the directors, thereby revealing detailed features on the striped domain structure, as shown in Fig. 4a. Figure 4b shows one representative normalized Raman map of the B_2g mode, measured in the (xx) polarization geometry, and Fig. 4c shows the corresponding AFM image. In addition, an R² map is also presented to assess the reliability of our fittings for constructing the director map in Fig. 4d, where there is little large deviation of R² from the ideal value +1. To convince ourselves, we also plot the fitting curves with the raw data of the polar plots at three representative points (Fig. 4e).

Fig. 4: Real-space visualization of the JT domains.

a A director map obtained by position-sensitive angle-resolved Raman spectroscopy. b Representative Raman map of the B_2g mode normalized by the A_g mode measured in the (xx) polarization geometry. c AFM image and d R² map over the same region. e Fitting curves for three representative points, denoted by red spots in a. All scale bars are 1 μm.

The widths of the stripe domains in Fig. 4a are clearly seen to be ~250 nm. The orientations of the directors are alternately aligned to the [100]_C and [010]_C axes. At the boundaries of the stripe domains, the magnitude of the director becomes vanishingly small, and we hardly detect an oblique orientation from the crystallographic axes. This is because the usual wall width of ferroelastic domains is less than tens of nanometers^47,48. According to the Abbe diffraction limit, our spatial resolution is limited by d = 0.61 × λ/N.A = 0.61 × 532 nm/0.9~360 nm (here, λ and N.A. indicate the wavelength of the light and numerical aperture). The small magnitudes of the directors result from the spatial averaging of two equal-portion domains due to a narrow width of ~150 nm, as estimated from the topographic images (see Fig. 4c and Supplementary Fig. 3).

High-energy regime

Until now we have discussed on the fundamental phonon modes associated with JT distortion and the corresponding domain structures. Meanwhile, the high-energy Raman modes (>800 cm⁻¹) were proposed as collective excitations of the OO, i.e., orbitons, and their presence was measured experimentally³⁵ stimulating academic discussion on more detailed identities and possible origins, such as mixed modes with both orbital and phonon characters due to the phonon-orbiton interaction^40,41, and multiphonon scattering activated by the Franck-Condon mechanism for a vibrationally excited self-trapped orbiton state^42,43. A consensus of the origin has not been established, but the orbital characteristic is included in the proposed mechanisms. Thus, the observation of the orbiton-related modes and their correlation with the underlying phonon modes provide useful information for understanding the domain structure of the local OO beyond the JT distortion.

With our analysis of the phonon modes in mind, we investigate the real-space distribution of the orbiton-related Raman peak intensities, which are involved in a local OO in our system. Figure 5a shows the zoomed-in Raman spectra for the high energy regime, normalized by the intensity of the fundamental A_g phonon mode. The inset figure is the Raman map for the B_2g mode measured in (yy) polarization geometry. Three representative positions were considered, as marked in the inset of Fig. 5a. Two orbiton-related peaks were clearly observed. The labeled A and B peaks show different tendencies in the intensity variation depending on the spatial position, suggesting the A and B peaks have A_g and B_2g symmetry, respectively. It is supported by the fact that the real-space intensity modulation of the orbiton-related B peak coincides with that of the fundamental phononic B_2g mode, as compared in Fig. 5b. Since the magnetic structure of LMO at low temperatures is strongly coupled to the orbital configuration via the superexchange interaction³⁸, the characterization of orbital-related domains offers a useful insight into the magnetic domain texture of LMO.

Fig. 5: Relation between the orientation of JT planes and orbiton-related peaks.

a Raman spectra in a higher-energy regime (> 800 cm⁻¹) at three representative positions. The peaks of labels B and A have B_2g and A_g symmetries, respectively. b Line profiles of the normalized Raman intensities of the fundamental B_2g mode (black) and the B peak of the orbiton (red). The dashed line in the inset indicates the corresponding positions. All data were obtained in the (yy) polarization configuration. (Insets) An identical Raman map for the normalized intensity of B_2g mode. The color scales vary from 0.8 to 1.5. The spatial scale bars stand for 1 μm.

Crystallographic structure

To elucidate the microscopic origin of the stripe domains, we characterized the crystal structure using x-ray θ-2θ scans and reciprocal space maps (RSMs). Figure 6a displays the θ-2θ scan exhibiting the (00 L)_PC diffraction peaks of the LMO film without noticeable impurity peaks (hereafter, ‘_PC’ stands for the pseudocubic indices of LMO film). From the θ-2θ scan, the out-of-plane lattice parameter of the LMO film was estimated to be 3.98(6) Å, which is almost the same as the bulk in-plane lattice parameters in terms of a pseudocubic unit cell³⁹. This implies that the orientation of the JT plane in our film is perpendicular to the substrate surface¹⁵.

Fig. 6: Structural analysis of the ferroelastic domains.

a A conventional θ-2θ scan along [00 L]_PC including the (001)_PC, (002)_PC, and (003)_PC LMO peaks. b ω-rocking curves measured along the [110]_C direction at the LMO (00 L)_PC reflections. c RSM for the (103)_PC reflections. d RSM for the (223)_PC reflections. The red arrows indicate peak shifts induced by the mosaic rotations. The black (red) dots in RSMs represent the diffraction peak positions of the M_c-type monoclinic unit cells without (with) mosaic rotation. e Schematic of a monoclinically distorted pseudocubic unit cell of LMO with (thick red) and without (pale red) mosaic rotation by δ. f Four possible variants of the monoclinic tilts. A solid (dashed) border indicates the top (bottom) surfaces of the unit cells. g A configuration of a⁺ and b⁻ unit cells with cooperative mosaic rotations and the top view. h The schematic of the a⁺/b⁻ domains with an W wall. i A configuration of the unit cells without mosaic rotation and the top view. j The schematic of the flat a⁺/b⁻ domains with an W’ wall. The black dashed lines in g–j represent domain walls. The cross-sectional schematics (bottom) in h, j correspond to line profiles along red solid arrows.

In Fig. 6b, the structural modulation is reflected in a ω-rocking curve measured along [110]_C at the LMO (001)_PC, (002)_PC, and (003)_PC reflections. In addition to strong central peaks, two satellite peaks were observed at the identical off-axis ω-angles of ±1.1(7) (°) regardless of the out-of-plane reciprocal positions, which suggests the mosaic rotation of the unit cell rather than structural periodicity along the [110]_C direction.

To clarify the crystal structure including in-plane lattice parameters, we obtained RSMs for two asymmetric reciprocal positions around the (103)_PC and (223)_PC peaks (Fig. 6c, d). To interpret the RSMs, we first considered a tilted pseudocubic unit cell (M_C-type monoclinic unit cell of which the c-axis tilts toward <100> _PC direction⁴⁹) of the orthorhombic LMO (Fig. 6e). Given a monoclinic angle β_PC with four-fold twins (denoted by a^± and b^±, respectively) imposed by the cubic symmetry of the STO substrate (Fig. 6f), the RSMs around the (103)_PC reflections result in three split peaks, corresponding to the black crosses in Fig. 6c. Two of them, corresponding to (103)_PC and ($\bar 1$03)_PC reflections, are shifted toward the –L or +L direction at an identical distance in the reciprocal space from the expected pseudotetragonal position estimated from the (00 L)_PC reflections, which is related to monoclinic angle β_PC. Additionally, the peak at the pseudotetragonal position corresponds to (013)_PC and (0$\bar 1$3)_PC reflections. The in-plane reciprocal positions give us information about the in-plane a_PC- and b_PC-axis lattice parameters.

Regarding the unit cell orientation, we further considered the mosaic rotation (δ) of the monoclinic unit cell with respect to the diagonal ([1$\bar 1$0]_C) direction, as shown in Fig. 6e, to interpret the remaining peaks (red crosses in Fig. 6c). The RSM in Fig. 6c reflects the mosaic tilts leading to the rotation of the peak positions, and further peak splits tracing along the polycrystalline arcs (the yellow lines in Fig. 6c). Because of the mosaic rotation, a clockwise (counterclockwise) rotation arises for the negative (positive) sign of the Miller indices, as clearly seen in the rotation of the (103)_PC and (0$\bar 1$3)_PC reflections. Subsequently, the overlapped peaks in the pseudotetragonal L position are split into (013)_PC and (0$\bar 1$3)_PC reflections by the mutually opposite rotations.

To ensure the model unit cell, we checked the RSM around the (223)_PC reflection. According to our models, two variant peaks should appear along the L-axis due to the monoclinic tilt, and two more peaks are created by mosaic rotations with different rotation directions along the polycrystalline arcs (Supplementary Fig. 4d). Figure 6d displays the resultant reflections around the (223)_PC reflection. The black crosses stem from the monoclinic unit cells with no mosaic rotation, and the red crosses correspond to the unit cells subjected to mosaic rotation. It is worth noting that the reciprocal positions of the (223)_PC and ($\bar 2\bar 2$3)_PC reflections are aligned to the pseudotetragonal L position. This arises from the fact that the two mosaic domains have the mirror symmetry of each other with respect to a {110}_PC plane, and should be continuously connected to minimize the elastic energy cost within the domain walls (see the schematics in Fig. 6g). In other words, the originally tilted c-axis vector is rotated so that the in-plane component is parallel to the axis of rotation for matching at the wall. Within the framework of our monoclinic models, we were able to determine the pseudocubic lattice parameters of the LMO film to be a_PC = 3.97(1) Å, b_PC = 3.84(1) Å, c_PC = 3.98(6) Å, β_PC = 88.2(1)°, and the mosaic rotation angles to be δ = 0.0(3)° for the upright unit cell and δ = 1.1(7)° for the unit cell subjected to the severe mosaic rotation.

Two types of domain and wall textures

Next, we link the structural information to the surface topography for the two monoclinic unit cells with and without mosaic rotation. To do that, we first try to interpret the orientations of domain walls based on the compatibility relation of ferroelastic domains⁵⁰,

$$\vec x^T[\mathop{\varepsilon}\limits^{\leftrightarrow} ({\mathrm{S}}_{\mathrm{1}}) -\mathop{\varepsilon}\limits^{\leftrightarrow} ({\mathrm{S}}_{\mathrm{2}})]\vec x = 0.$$

(9)

Here, we consider two spontaneous strain tensors $\mathop{\varepsilon}\limits^{\leftrightarrow} ({\mathrm{S}}_{\mathrm{1}})$ and $\mathop{\varepsilon}\limits^{\leftrightarrow} ({\mathrm{S}}_2)$ where strain states S₁ and S₂ come from an identical prototype in high symmetry. And, $\vec x$ is the position vector of the domain wall. According to the quadratic Eq. (9), there are two permissible domain walls between two distinct strain states (e.g., a⁺ and b⁻; see Fig. 6g–j, Supplementary Fig. 5 and Note); one is the vertical mirror plane along [1$\bar 1$0]_C called W wall and the other is the inclined plane along [110]_C called W’ wall. These domain wall orientations in M_C-type pseudocubic unit cells are contrasted with domain and wall configurations along <100>_PC orientations in rhombohedral PbZr_1-xTi_xO₃ and manganites^51,52,53 because the M_A-type monoclinic pseudocubic unit cells (of which the c-axis tilts toward <110>_PC direction) result in domain walls along <100> _PC by azimuthal rotation of 45° relative to the present M_C-type case. The W wall is formed when the monoclinic tilt directions collide at the wall. To minimize the interfacial elastic energy at the W wall, a⁺ and b⁻ domains are subjected to a clockwise or counterclockwise mosaic rotation (δ = ±1.1(7)°) around the [1$\bar 1$0]_C axis in a symmetrical manner (Fig. 6g, h). For a perfect match at the wall, the mosaic rotation angle should be the monoclinic angle (90 − β_PC) divided by $\sqrt 2$, that is ~1.2°. However, this mosaic rotation inevitably makes the interfacial energy between the film and the substrate more unstable. So there is a trade-off problem, and the trade-off will depend on the film thickness. The thickness of the LMO film studied here is ~150 nm, but we cannot observe the protruding surface topographic feature when the film is as thin as 35 nm (Supplementary Fig. 6).

Meanwhile, the W’ wall is formed when the monoclinic tilt direction continuously changes across the domain wall into a head-to-tail type, as shown in Fig. 6i, j. It doesn’t require a mosaic rotation (δ = 0.0(3)°) to relieve the interfacial energy. These expected properties of ferroelastic domain walls show a good agreement with the observed topographic features. We conclude that the regions with the bumpy surface in the AFM images consist of a^±/b^± domains with W walls, reflecting the angle of the slope of the surface, ±1.(0)° (Supplementary Fig. 3), similar to the value of mosaic rotation measured by x-ray diffraction. Meanwhile, the flat regions have W’ walls (Fig. 6j). It is noteworthy that the W’ wall is inclined to the substrate by the amount of 56° considering the strain state and the compatibility relation, and the W wall is perpendicular to the substrate (Fig. 6h, j and Supplementary Fig. 5).

The periodicities of stripe domains obtained by AFM and Raman spectromicroscopy are quite different, and they are definitely not one-to-one correspondence (Fig. 1e, g). One is topographically flat and the other has an uneven surface. The flat area (represented by the solid circle in Fig. 1) shows the periodic stripe domains on the Raman map in B_2g mode, while the topographically modulated area (denoted by the dashed circle in Fig. 1) has no distinctive stripes. We estimated that the domain widths of the surface modulated region to be ~150 nm (300 nm in periodicity) based on AFM line profile data shown in Supplementary Fig. 3, while the domain widths of the flat region to be ~250 nm (~500 nm in periodicity) based on the line profile data from the Raman map shown in Fig. 5b. The seeming inconsistency between the surface topography and the Raman map can be understood based on the spatial resolution of our Raman imaging. The two possible domain textures with different mosaic orientations have the different typical widths of domains. One happens to be larger and the other is smaller than the spatial resolution of our Raman imaging. According to the Kittel’s law, the typical size of domains is determined by the competition between the domain energy and the domain wall energy, and the interfacial energy between substrate and film should be also considered. Although it is natural that the typical domain widths for two domain/wall textures are different, further study is required for quantitative theoretical estimates.

Redistribution of the JT domains by mechanical writing

Based on the ferroelastically-formed domain structure, we demonstrate that the distribution of the JT domains could be manipulated by applying a mechanical force, using a AFM tip. We exerted a strong frictional force associated with a normal force of ~150 µN corresponding to ~100 GPa (approximating the tip-surface contact to be a disk with 20 nm in radius). We scanned along two different directions, [100]_C and [1$\bar 1$0]_C, to test whether the trailing effect could align the domains and walls (Fig. 7a).

Fig. 7: Mechanical writing of JT domains by AFM tip.

a Area scans by a tip imposed by a normal force of ~150 μN. The black lines indicate the fast scan direction of a tip along 0° and −45° with respect to [100]_C. AFM images of the LMO film before (b) and after (c) the mechanical writing. The black boxed image is a zoomed-in image of the region marked in c. Raman maps of the B_2g mode before (d) and after (e) the writing over the same region. Transparent dashed boxes indicate the scanned regions. The black arrows indicate representative lines where the bundles of the W-wall-type domains collide other orthogonal bundles. Scale bar represents 1 μm.

Figure 7 shows the effect of mechanical writing on the topography and Raman signals of the B_2g mode. The mechanical scans erased the W’-wall-type domains (δ = 0°), but the regions became filled with W-wall-type domains (δ = 1.1°). It is attributed to the strong frictional force of the tip, which induces the mosaic rotation and/or forces the monoclinic tilt direction. We note that the written domain walls become denser after the writing. The scan direction also affects the density; the diagonal scan induces a higher density of domain walls with a lower height modulation. The directional control of domain walls doesn’t seem to be fully deterministic. The nucleated walls tend to propagate in a straight line, but some have to break the continuation due to the intervention of adjacent configurations. Similar mechanical writing has been accomplished in ferroelastic WO₃ films⁴⁸, and the mechanical load used in the current study is much larger than before. We cannot exclude the possibility that the applied stress would only affect the top surface of the film and thus the depth of the writing would be shallower than the film thickness. Since the scan along [100]_C is parallel to the monoclinic tilt direction, the writing effect could reach deeper in the film and the domain width might be wider than the [1$\bar 1$0]_C case.

It is interesting to note that a distinctively enhanced Raman signal arises along the boundaries where the bundles of the W-wall-type domains collide with another orthogonal bundles (indicated by the arrows in Fig. 7c, e). These supermacro domain boundaries are expected to undergo severe local lattice deformation with a large strain and strain gradient and might also be involved in orbital frustration. We note that a strain gradient breaks local inversion symmetry, thereby flexoelectric polarization can be induced⁴⁸. It has been reported that Raman signals are enhanced at ferroelectric domain walls⁵⁴ as a consequence of local strain and/or electric fields that can improve Raman scattering efficiency.

We successfully observed the real-space distribution of JT domains in an epitaxial LMO thin film using confocal Raman spectromicroscopy. Based on the symmetry analysis of the Raman modes, we could identify and visualize the local populations and orientation of the JT domains, and the concomitant OO distribution was identified. A detailed structural analysis revealed the formation of ferroelastic domains, which becomes a framework for the JT domain walls under the consideration of the ferroelastic compatibility relation. We also performed mechanical manipulation of the JT domains. Our findings offer a promising approach for the observation of orbital domains, and provide useful information for understanding JT-mediated physics.

Methods

Synthesis of the LMO thin films

LMO thin films were deposited on a (001)_C-oriented Nb-doped (0.5 wt%) STO substrate (CrysTec GmbH) using pulsed laser deposition with a KrF excimer laser (λ = 248 nm). Growth was made at a temperature of ~900 °C in the Ar and O₂ gas mixture of 0.1 Torr with a ratio of 9:1. Laser fluence was ~ 0.5 J/cm² with a repetition rate of 10 Hz. After deposition, the films were cooled down to room temperature at a rate of 10 °C min⁻¹ in vacuum environments (~10⁻⁶ Torr) to resolve the cation vacancy issue, by reducing excess oxygen.

X-ray diffraction analysis of the LMO thin films

To investigate the crystal structure of the LMO film on Nb-doped STO, x-ray θ-2θ scans and RSMs were measured at room temperature using a high-resolution x-ray diffractometer (PANalytical X’pert-PRO MRD) with Cu Kα₁ radiation (λ = 1.5406 Å).

Characterization of surface morphology and mechanical writing

The surface topographies of the LMO thin film were characterized using AFM (MFP-3D Infinity, Oxford Instruments, UK) with Si tips (HQ:NSC35, MikroMasch) with a spring constant of 5–20 N/m at scan rates of 5–10 μm/s. The force load during the surface measurement was less than ~80 nN. For mechanical writing, we used Sb-doped Si tips (RTESPA-525, Bruker) with a spring constant of 200 N/m.

Raman spectroscopy

We examined the phonon modes of the LMO thin films using the confocal Raman spectroscope alpha 300R (WITec GmbH, Germany). Raman spectra were excited using a Nd:YAG 532 nm laser, and were recorded in backscattering geometry using a UHTS 600 spectrometer. A laser power density of ~20 kW/mm² at the sample surface was used. We used an integration time of 0.25 s at each position to obtain the Raman maps of the fundamental phonon modes. To get line profiles of the high-energy regime, an integration time of 30 s was required at each position.