Search tips
Search criteria 


Logo of ncommsLink to Publisher's site
Nat Commun. 2017; 8: 615.
Published online 2017 September 20. doi:  10.1038/s41467-017-00710-5
PMCID: PMC5607007

Controlled manipulation of oxygen vacancies using nanoscale flexoelectricity


Oxygen vacancies, especially their distribution, are directly coupled to the electromagnetic properties of oxides and related emergent functionalities that have implications for device applications. Here using a homoepitaxial strontium titanate thin film, we demonstrate a controlled manipulation of the oxygen vacancy distribution using the mechanical force from a scanning probe microscope tip. By combining Kelvin probe force microscopy imaging and phase-field simulations, we show that oxygen vacancies can move under a stress-gradient-induced depolarisation field. When tailored, this nanoscale flexoelectric effect enables a controlled spatial modulation. In motion, the scanning probe tip thereby deterministically reconfigures the spatial distribution of vacancies. The ability to locally manipulate oxygen vacancies on-demand provides a tool for the exploration of mesoscale quantum phenomena and engineering multifunctional oxide devices.


Oxygen vacancies (Vo) are elemental point defects in oxides, and they generally function as mobile electron donors. At sufficiently high concentrations, they can disturb the ground state and promote emergent functional phenomena, such as superconductivity1, ferromagnetism2, 3, metal-to-insulator transition4 and interface conductivity5. In addition, the distribution and dynamics of Vo play a central role in many oxide-based energy and memory applications6, 7. Hence, the ability to manipulate Vo provides an opportunity to study the evolution of emergent phenomena and control numerous functionalities, which are essential for developing next generation oxide devices8, 9.

Traditionally, modification of the vacancy concentration has been carried out by annealing oxides at high temperature in reducing/oxidising atmosphere. Recent works have shown that at room temperature, the vacancy concentration can instead be locally changed by employing an electrical bias from a scanning probe microscope (SPM) tip4, 10. Such local manipulation enables reversible nanoscale control of metal-insulator transitions4 and the modulation of interface conductivity11. However, in addition to the change in the vacancy concentration, the application of bias through an SPM tip is often accompanied by charge injection12, 13 and the formation of protons/hydroxyls14, which complicates the practical implementation of this approach.

Recently, it has been shown that the mechanical force from the SPM tip can also deplete Vo from the contact region15, 16. However, because this mechanical depletion is less understood, the feasibility of using the force from an SPM tip as an active means of manipulating Vo has not been explored. An explanation for the mechanical depletion of Vo has been proposed through the so-called piezochemical coupling mechanism1517, which involves the converse Vegard effect and the flexoelectric effect. The converse Vegard effect accounts for the decrease in vacancy concentration due to force-induced lattice compression. Meanwhile, the flexoelectric effect considers the electromigration of positively charged Vo under the stress-gradient-induced flexoelectric field. Notably, by definition, the flexoelectric effect refers to the generation of an electrical polarisation by stress-gradients18, 19. In this context, Yudin and Taganstev suggested that the flexoelectric field is not a macroscopic electric field19; instead, it should be understood as a pseudo-internal field that polarises a medium in the presence of a stress-gradient2023. Hence, it remains unclear how this flexoelectric field acts on Vo.

To be able to mechanically manipulate Vo, we must, therefore, thoroughly understand how they respond to the mechanical force from the SPM tip. Subsequently, we have to devise a strategy to employ this force to increase as well as decrease the vacancy concentration in a controlled manner.

Here, we demonstrate a controlled manipulation of Vo using the mechanical force from an SPM tip. As a model system, we used a homoepitaxial thin film of SrTiO3, which is an archetypal quantum paraelectric oxide with well-known flexoelectric coefficients24, 25. Using a Kelvin probe force microscopy (KPFM)-based imaging scheme, we show that besides pushing Vo away in the vertical direction, the force promotes the lateral transport of vacancies during the movement of the tip. Using phase-field simulations, we argue that this mechanical redistribution of Vo is driven by the depolarisation field associated with the stress-gradient-induced flexoelectric polarisation. The depolarisation field underneath the tip promotes the vertical migration of vacancies, whereas around the contact edge, it traps Vo that can move with the tip. Furthermore, we demonstrate that altering the tip geometry can tailor this depolarisation field to preferentially allow the lateral transport of Vo. This approach enables a controlled spatial modulation of vacancies.


Probing oxygen vacancies with Kelvin probe force microscopy

We start by discussing the concept of probing distribution of Vo with the KPFM technique (Fig. 1a–c). KPFM is a surface sensitive technique that measures the contact potential difference (CPD) between the tip and SrTiO3 (STO) surface. The CPD is defined as the offset in the respective vacuum energy levels (in units of V)26. However, as schematically illustrated in Fig. 1b, this CPD would change if the STO surface contains Vo, which can be locally accumulated by an electrical poling. This change in the CPD can originate from the workfunction of STO27, chemical dipoles (including Vo-electron pairs), and their orientation28. Moreover, the tip bias used during the KPFM measurement has been argued to influence the measured CPD29. The contribution of each factor cannot be individually separated, which inhibits a quantitative determination of the vacancy concentration with the KPFM technique. Nonetheless, as shown in Fig. 1c, the KPFM contrast across this Vo-rich region can be argued to scale proportionally with the vacancy concentration27. To establish a proof of concept, in the following we elaborate on the application of KPFM technique to study the diffusion of Vo using a 14 and 120-unit cell (uc)-thick STO film (details of the film are provided in the Methods section and in Supplementary Fig. 1).

Fig. 1
Studying diffusion characteristics of oxygen vacancies with KPFM. a Sketch of the sample geometry and KPFM measurement architecture. b Schematic illustration of the contact potential difference (CPD) contrast between the pristine (referred to as 1) and ...

Figure 1d, e show KPFM images after poling the pristine surface of STO films with a tip bias of −5 V. The images of the 14-uc (120-uc)-thick STO film were taken 20 (16) minutes and 380 (360) minutes after poling. An image contrast is clearly visible across the poled region, which is appearing as a dark rectangular patch. This implies surface charging, either by injected charges, protons/hydroxyls, and/or Vo 30, 31. The injected charges and protons/hydroxyls are expected to decay with a timescale, which is independent of the STO thickness31. In contrast, because of diffusion or surface reaction that enables the recombination of Vo with oxygen from the ambient, the vacancy concentration is expected to decrease with a pronounced thickness-dependent timescale32. Figure 1d, e indeed show that the contrast across the poled region diminishes over time, and the decrease is largest in the 14-uc-thick STO film. Therefore, we argue that the surface charging is caused by Vo, which either undergo diffusion or recombine with oxygen from the ambient.

We can distinguish the dominating mechanism that causes the vacancy concentration to decay over time by analysing the time-dependency of the degree of equilibrium, S(t) that can be calculated from the KPFM images (see Supplementary Note 1). The time-dependency of S(t) describes how the surface equilibrates after the electrical poling. S(t) will exhibit a semi-parabolic (linear) time-dependency if diffusion (surface reaction) is the dominating mechanism27, 32. Figure 1f plots S(t) as a function of time. Evidently, the time-dependencies of S(t) are semi-parabolic for both films–implying that the diffusion of Vo causes the vacancy concentration to decrease over time.

Arguably, during the poling, the applied tip bias perturbs the equilibrium Vo-distribution of the entire film including the surface, which equilibrates through the diffusion of Vo. This diffusion occurs both along the out-of-plane and in-plane directions. However, due to the high surface sensitivity of the KPFM technique, the surface-bulk diffusion of Vo predominantly affects the time evolution of the KPFM contrast, and thus the time-dependency of S(t). Notably, during this surface-bulk diffusion, the repulsive vacancy–vacancy interaction inhibit Vo to migrate independently along the out-of-plane direction33. Thus, the time-dependency of S(t) effectively describes how the perturbed volume under the poled area (in Fig. 1d, e) equilibrates (see Supplementary Note 1 for a detailed discussion). Naturally, this volume would be smaller in the thinner STO film, and thus would equilibrate faster. This explains why the KPFM contrast (S(t)) diminishes (grows) more rapidly in the 14-uc thick STO than in the 120-uc thick STO film (Fig. 1d–f).

Following the rationale above, we fit the time evolution of S(t) with Fick’s 2nd law of diffusion (solid lines in Fig. 1f). From this fitting, we obtained the diffusion coefficient D = 9.4(3) × 10−19 and 3.8(2) × 10−18 cm2 s−1 for the 14-uc and 120-uc-thick STO film, respectively. These values are well in the range of the bulk value D bulk  10−(17±3)  cm2 s−1 (300 K-extrapolated)34, which validates the conceptual schematic depicted in Fig. 1c. Subsequently, we utilised this correlation between the vacancy concentration and KPFM signal to image the vacancy redistribution under an applied bias and force.

Oxygen vacancy redistribution under applied bias and force

In this section using the 120-uc-thick STO film, we compare the response of Vo to the electrical bias and force from the SPM tip. For this study, we employed a sharp tip with a radius of curvature = 25 nm. As depicted in the KPFM image (Fig. 2a), we scanned the left side of a formerly Vo-enriched region with an increasingly positive tip bias (0–5 V) at a contact force = 0.6 µN. Meanwhile, the right side was scanned with an increasing contact force (0.6–8.5 µN) by holding the tip at the ground potential. In either case, the tip crossed the border between the formerly Vo-enriched and pristine regions. The central part of the KPFM image is uniformly dark, which indicates a uniform vacancy distribution. However, a change in the image contrast, which stems from the Vo-redistribution is visible within both electrically and mechanically scanned areas. Notably, these scans together with acquiring the KPFM image took about 20 min. Because of their ultra-slow decay as evident from Fig. 1e, during this time Vo can be assumed to be kinetically frozen. This implies that the Vo-redistribution is due to the applied bias and force.

Fig. 2
Characterisation of oxygen vacancy redistribution by applied bias and force. a The KPFM image after electrical and mechanical scans were performed across borders between the Vo-enriched and pristine regions. The schemes for electrical ...

To quantify this Vo-redistribution, we constructed the vacancy concentration map from Fig. 2a by defining the normalised vacancy concentration, NVC=Vb-VVb-Vmin. Here, V b and V min represent the baseline and minimum KPFM signals, which were extracted from areas indicated by the black (top-right corner) and white (centre) squares, respectively. The reconstructed NVC map is shown in Fig. 2b. The positive (negative) NVC refers to a higher (lower) vacacny concentraion compared to the intrinsic surface concentration of Vo (NVC = 0, at the top-right corner).

Profiling the NVC map along lines E and M indicates that both the applied bias and force depleted the formerly Vo-enriched region (Fig. 2c). We find that this depletion can be modelled with a phenomenological Boltzmann sigmoid function in the following form,


where A 1, A 2, X o, and M refer to the initial value, final value, center, and decay rate, respectively. The solid lines in Fig. 2c are fits to this function. Table 1 lists the best-fit parameters. This functional analysis suggests that the applied positive bias and force have a same qualitative effect along the lines E and M, respectively. By comparing the X o ( = 0.5A 1) values from the fits, we obtain a force-voltage equivalence factor of 0.4 V µN−1. The Boltzmann sigmoid fits in Fig. 2c also exhibit an onset. Defining this onset as 0.9A 1 (indicated by stars), on an ad hoc basis, yields a threshold voltage (V th) and force (F th) of 1.3 V and 2.7 µN or equivalently 1.1 V, respectively. These threshold values are consistent with the activation barrier potential ( = 1−1.4 V) for the electromigration of Vo in STO35, 36, which further validates the force-voltage equivalence relationship.

Table 1
Best-fit parameters from the Boltzmann sigmoid fit to the normalised vacancy concentration (NVC) data in Fig. 2c

Surprisingly, electrical and mechanical scans yield the opposite image contrast across borders in the NVC map. As shown in Fig. 2d, along line EL, NVC monotonically decreases with the increasing bias. However, along line ML, NVC varies non-monotonously: first, it increases as the force increases and then gradually drops; the horizontal line in Fig. 2d marks the background (NVC = 0.1 at 0 V and 0.6 µN). We, therefore, conclude that the electrical scan depletes the pristine surface outside the left border, while the mechanical scan enriches the pristine surface outside the right border.

The force-induced vacancy-enrichment of the pristine region implies that some of the depleted Vo moved laterally with the tip across the border during the mechanical scan. Additional experiments elaborate that the mechanical scan does not alter the background within the pristine region, and thus rule out the formation of Vo and triboelectric charging of the STO surface during the mechanical scan (see Supplementary Note 2). To quantify the lateral motion of Vo, in Fig. 2e we, therefore, show the NVC with the background (NVC at 0.6 µN) subtracted (i.e., ΔNVC) along lines M and ML. Because KPFM is a surface-sensitive technique, the remnant Vo on the surface contribute to the NVC map after Vo-redistribution. Also, the possibility of recombination with oxygen from the ambient during the lateral motion of Vo can be ignored based on the following considerations. The influence of the surface reaction process, which could facilitate this recombination is negligible in our film. Furthermore, the use of a grounded tip during the mechanical scan rules out the bias-induced amplification of this surface reaction process37. Hence, for a particular force, the ratio between the net gain (along ML) and the net drop in NVC (along M) represents the fraction of depleted Vo that moved laterally with the tip. By comparing the maximum net gain ( = 0.1), and corresponding drop ( = 0.3) at force = 5.7 µN, we estimated that only approximately 1/3 of the depleted Vo moved with the tip along the surface, while the rest migrate into the bulk.

Overall, our main observations are as follows. First, the contact force is bifunctional: it depletes the formerly Vo-enriched region and simultaneously enriches the pristine surface. Second, the mechanical Vo-redistribution involves a predominant surface-bulk migration and a relatively weaker lateral motion of Vo with the tip.

Modelling the mechanical redistribution of oxygen vacancies

To understand the mechanical redistribution of Vo we first considered two mechanisms–the converse Vegard effect and the flexoelectric effect15, 16. Recently, the magnitude of these two effects under an applied force from SPM tip has been compared in a PbTiO3 (PTO) thin film38. This study suggests that the converse Vegard effect is much weaker than the flexoelectric effect. Notably, both the PTO and STO have comparable flexoelectric and Vegard coefficients25, 38, 39, which determine the relative contributions of these two effects for a given force. Based on these considerations, we thus conclude that the flexoelectric effect predominantly causes the mechanical redistribution of Vo, and the contribution from the converse Vegard effect is marginal.

For gaining a mechanistic understanding of this flexoelectric effect-driven Vo-redistribution, we performed phase-field simulations; whereby we incorporated the flexoelectric effect40 and coupled the time-dependent Ginzburg-Landau and the Nernst-Planck equations41, 42. The simulation was performed assuming that the STO is paraelectric (see Supplementary Note 3). Unlike in the experiment, the SPM tip was assumed to be static, and following the Hertzian model the contact radius was calculated to be 8 nm for a contact force of 4 µN. Initially, Vo were assumed to be homogeneously distributed over the entire STO thickness, instead of being localised on the surface. Despite these oversimplified assumptions, our simulation still provides a qualitative insight into how Vo respond to mechanical stimuli. A detailed explanation of our model and a discussion on the redistribution of Vo under a positive tip bias are included in Supplementary Notes 4, 5.

The stress-gradient from the SPM tip locally polarises STO through the flexoelectric effect24, as evident from the out-of-plane polarisation vector map in Fig. 3a. Since in the simulation the STO film is assumed to be in the paraelectric phase, the polarisation in Fig. 3a purely stems from the flexoelectric effect. This flexoelectric polarisation reaches a maximum ( ~ 0.04 C m−2) underneath the tip and spatially varies both in magnitude and direction. The resulting polarisation bound charge and the associated depolarisation field accordingly redistribute vacancies around the tip-STO contact region. This redistribution can be visualised from the simulated NVC map in Fig. 3b, which plots the in-plane distribution of Vo at the surface. Clearly, the vacancy concentration is decreased (enhanced) underneath the tip (around the contact edge).

Fig. 3
Phase-field simulations. a Simulated out-of-plane (z-x plane) vector map of flexoelectric polarisation induced by the scanning probe microscope (SPM) tip under a static contact force of 4 µN. Arrows denote the direction of the polarisation ...

An intuitive understanding of the simulated Vo-redistribution can be gained from the component-resolved distribution of the depolarisation field. As shown in Fig. 3c, the z-component, Ezdep, points downward (upward) below the tip (around the contact edge), whereas the x-component, Exdep, exhibits a parentheses-like structure: a node at the contact point and antinodes around the contact edge (Fig. 3d). Effectively, Exdep points inwards, as indicated by the white arrows in Fig. 3d. The y-component (not shown), Eydep, forms an analogous structure to Exdep but is rotated by 90° in the x-y plane. The decrease in the vacancy concentration underneath the tip can be attributed to the downward Ezdep component, which moves positively charged Vo from the surface into the bulk. In contrast, the combination of the upward Ezdep and inwardly directed Ex,ydep components favours the accumulation of Vo around the contact edge, yielding an increase in the vacancy concentration.

The simulation results qualitatively explain the characteristics of the mechanical Vo-redistribution in Fig. 2. During a mechanical scan, a spatially extended and strong downward Ezdep field acts over a larger fraction of Vo underneath the tip, which results in a dominant surface-bulk migration. In contrast, a small fraction of Vo becomes effectively trapped within a shallow annular region around the contact edge by the upward Ezdep and inward Ex,ydep. These trapped vacancies can move laterally with the tip from the Vo-enriched region to the pristine region during the tip’s lateral motion. Therefore, the contrasting roles of the depolarisation field underneath the tip and around the contact edge corroborate both the dominant surface-bulk migration and the relatively weaker lateral motion of Vo with the tip along the surface.

Notably, while scanning, the tip redistributes Vo regardless whether it moves from left-right (trace) or right-left (retrace). Thus, the mechanical redistribution of Vo should be understood as an average response of Vo to the force applied during the trace and retrace. Moreover, the scanning velocity would influence the Vo-redistribution–longer the tip spends in contact with STO, the larger number of Vo it would redistribute. To check whether these factors could contribute to the weaker lateral motion of Vo, we performed additional experiments, whereby we applied force only during the trace, and varied the scanning velocity (see Supplementary Figs. 911). These experiments also yielded a weaker lateral motion but a stronger surface-bulk migration of Vo–highlighting the dominating influence of depolarisation field underneath the tip. In the following, we illustrate that the depolarisation field around the tip-STO contact junction can be tailored in favour of the lateral motion of Vo, which enables controllably manipulating the vacancy distribution.

Controlled spatial modulation of oxygen vacancies

The basic concept of tailoring the SPM tip-induced depolarisation field can be understood with the aid of Fig. 4a, which compares the simulated surface deformation profiles under a static load of 4 µN using two different tip geometries. The upper panel corresponds to the spherical tip (contact radius = 8 nm) that is used in Fig. 3, and the lower panel corresponds to a flat-ended tip (contact radius = 15 nm). Compared to the spherical one, the flat-ended tip usually imparts a weaker stress on the STO surface underneath. Since the downward Ezdepscales with the stress-gradient underneath the tip, the flat-ended tip induces a very small downward Ezdep (Fig. 4b, upper panel). In contrast, the lateral deformation (indicated by curved arrows in Fig. 4a), which controls the depolarisation field around the contact edge, is alike for both geometries. Consequently, the flat-ended tip induces a depolarisation field distribution around the contact edge (Fig. 4b) similar to that of the spherical tip in Fig. 3c, d. Additionally, an enhanced contact radius enlarges its spatial extent. A selective suppression of Ezdep underneath a flat-ended like tip should significantly reduce the surface-bulk Vo migration. Meanwhile, the extended depolarisation field around the contact edge should improve the lateral transport of Vo.

Fig. 4
Controlled manipulation of oxygen vacancies. a Simulated surface deformation profiles under a spherical (upper panel) and flat-ended (lower panel) tip for a static contact force of 4 µN. b Simulated in-plane distribution of the z-component, ...

To validate our proposition, we performed experiments with a sharp and blunt tip. The estimated radius of curvature of this blunt tip is larger than 200 nm (see Supplementary Fig. 7). Thus, it effectively yields a flat contact junction underneath the tip. We used the 120-uc thick STO film in these experiments. Figure 4c, d show the NVC maps after mechanical scans were performed with these tips at a contact force of 9.5 µN. Notably, we scanned both the left and right boundaries between the Vo-enriched and pristine regions.

To compare the sharp and blunt tip-induced Vo-redistribution we profiled the NVC maps, as indicated by vertical lines in Fig. 4c, d. Figure 4e, f show the corresponding NVC profiles. The overlapping NVC profiles along lines M1/M4 and M2/M3 demonstrate that the Vo-redistribution is reproducible for both tip geometries. However, the response of Vo to the applied force from these two tips are clearly different. The NVC profiles in Fig. 4e exhibit a maximum drop in NVC of Δmaxdec=-0.75 along lines M2 and M3 but no appreciable increase in NVC (Δmaxinc) along lines M1 and M4. This implies that the sharp tip strongly depletes the Vo-enriched regions but barely enriches the pristine regions. This result is in qualitative agreement with that in Fig. 2e, which shows that the fraction of the depleted Vo that laterally move with the tip progressively decreases for applied forces larger than 6 µN. The NVC profiles in Fig. 4f, however, exhibit a maximum drop in NVC of Δmaxdec= −0.25 along lines M2 and M3 and increase by Δmaxinc= + 0.2 along lines M1 and M4. This implies that approximately 80% of the depleted Vo laterally moved with the blunt tip. The strong reduction (improvement) of Δmaxdec (Δmaxinc) thus confirms an active suppression of the out-of-plane migration and a simultaneous enhancement in the lateral transportation of Vo, during scans with the blunt tip. Overall, the ability to deterministically move Vo constitutes the first experimental demonstration of a controlled manipulation of Vo in an oxide and the resulting two-dimensional spatial modulation.


Through a combined experimental and theoretical approach, we demonstrated the flexoelectricity-mediated controlled manipulation of oxygen vacancies by the mechanical force from an SPM tip. A deterministic reconfiguration of spatial vacancy profile provides control over the electron density and related electronic correlation effects. This could enable, using an SPM-based all-in-one platform, the investigation of mesoscale quantum phenomena in oxides43. Ultimately, this creates the opportunity for developing mechanically sketched oxide devices, and ambipolar mechanical control of device functionalities such as electroresistance states. The voltage-free operation of the SPM tip would thereby eliminate the possibility of surface charging. At this point, we want to emphasise that flexoelectricity is a universal phenomenon, which can occur in any dielectric18, 44. However, the flexoelectric coefficients of few oxides, such as SrTiO3 and BaTiO3, are currently known45. Thus, our work should motivate the study of flexoelectricity in other oxides.

Broadly speaking, our KPFM-based imaging approach offers a time-efficient way of characterising the activation barrier potential for oxygen vacancy migration at room temperature to complement the conventional Arrhenius analysis46. Combined with the feasibility of determining the diffusion coefficient, this technique could thus become an essential metrology tool for oxide-based energy and memory research. Furthermore, our theoretical model that couples the phase-field simulations to the Nernst-Planck equation, can be employed to elucidate how depolarisation fields cause oxygen vacancies., electrons, and holes to redistribute. Therefore, the model can be extended to the study of emergent problems such as the domain wall conductivity, high electrical conductivity of morphotropic phase boundaries, and leakage current in ferroelectric oxides4750.


Thin film growth

SrTiO3 thin films were homoepitaxially grown on TiO2-terminated Nb:SrTiO3 (0.5% wt. doped) substrates using pulsed laser deposition technique. The growth dynamics and thickness were monitored by in-situ reflection high energy electron diffraction technique. The depositions were performed at 1000 °C and using an oxygen partial pressure of 5 × 10−7 torr. After deposition, films were annealed at 800 °C for an hour in a 1 torr oxygen atmosphere and subsequently cooled down to room temperature at a cooling rate of 20 °C per minute.

Kelvin probe force microscopy

KPFM measurements were carried out using the Asylum Research Cypher SPM at room temperature and under ambient conditions. Pt/Ir-coated metallic tips (NANOSENSORS™ PPP-NCHPt) with a nominal spring constant ≈ 40 N/m were used for electrical/mechanical scans and KPFM imaging. The KPFM measurements were obtained in the non-contact mode using a lift height of 30 nm and the typical scan parameters used are as follows: V ac = 1 V (peak-to-peak), f resonance = 250 kHz, and scan rate = 1 Hz. Before each experiment, the spring constant of the cantilever was accurately determined from force–distance measurements and thermal tuning methods. The contact force during mechanical scans was varied accordingly by controlling the set-point voltage.

Theoretical modelling

To model the oxygen vacancy redistribution by mechanical force, we performed phase-field simulations by coupling the time-dependent Ginzburg-Landau and Nernst-Planck equations.



In Eq. (2), P is the polarisation vector, L is the kinetic coefficient and F is the total free energy of the system, which includes Landau, electric, gradient and flexoelectric energy contributions. In Eq. (3),Vo, DVo, μVo, and ϕ denote the concentration, diffusion coefficient, and mobility of the oxygen vacancies and the electric potential, respectively. Detailed descriptions of the energy functional F and other relevant simulation parameters are presented in the Supplementary Note 4.

To solve Eqs. (2) and (3), the system was discretised into 100Δx × 50Δy × 500Δz grid points to implement the semi-implicit Fourier method45. The 120-uc-thick STO film is simulated to be 480Δz in thickness, while the substrate and air are each 10Δz in thickness. The parameters for the total free energy of STO were adopted from work of Y.L. Li et al.46. Following the work of R. Moos et al.47, the initial concentration of oxygen vacancies in STO is assumed to be 3.66 × 1014 cm−3 based on our thin film growth conditions. In addition, the diffusion coefficient is assumed to be constant with regard to pressure and calculated as 1.23 × 10−15 cm2 s−1 at room temperature47. A sufficiently long simulation time is used to ensure that the induced polarisation and oxygen vacancy concentration reach a quasi-steady state. At each time step, the electrostatic and elastostatic equilibrium equations are solved under the electric short-circuit48 and mechanical mixed boundary conditions49, respectively.

Data availability

The data that support the findings of this study are available from the corresponding authors upon reasonable request.

Electronic supplementary material


This work was supported by the Institute for Basic Science in Korea (Grant No. IBS-R009-D1). B.W. and L.-Q.C acknowledges the support by the National Science Foundation under the grant number DMR-1410714 and by the Penn State MRSEC, Center for Nanoscale Science, under the award NSF DMR-1420620. Y.C. and S.V.K were supported by the U.S. DOE, Office of Basic Energy Sciences (BES), Materials Sciences and Engineering Division (MSED) under FWP Grant No. ERKCZ07 (Y.C., S.V.K.). A portion of this research was conducted at the Center for Nanophase Materials Sciences, which is a DOE Office of Science User Facility. We would like to thank Prof. Jong-Gul Yoon and Prof. Jin-Seok Chung for discussions. We also acknowledge Dr Luke Sandilands and John Henry Gruenewald for carefully proofreading the manuscript.

Author contributions

Author contributions

S.D. conceived and planned this project under the direction of T.W.N. S.D. grew the films and performed structural characterisation with assistance from M.K. S.D. carried out KPFM measurements assisted by Y.J.S., S.M.Y., and L.F.W. and performed data analysis. B.W. and Y.C. performed theoretical modelling under the direction of L.-Q.C. and S.V.K. M.R.C. performed SEM measurements. S.D. and T.W.N. wrote the manuscript with inputs from all authors.


Competing interests

The authors declare no competing financial interests.


Electronic supplementary material

Supplementary Information accompanies this paper at doi:10.1038/s41467-017-00710-5.

Publisher's note: Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Contributor Information

Saikat Das, moc.liamg@78sad.takias.

Tae Won Noh,


1. Kang HJ, et al. Microscopic annealing process and its impact on superconductivity in T′-structure electron-doped copper oxides. Nat. Mater. 2007;6:224–229. doi: 10.1038/nmat1847. [PubMed] [Cross Ref]
2. Rice WD, et al. Persistent optically induced magnetism in oxygen-deficient strontium titanate. Nat. Mater. 2014;13:481–487. doi: 10.1038/nmat3914. [PubMed] [Cross Ref]
3. Santander-Syro AF, et al. Giant spin splitting of the two-dimensional electron gas at the surface of SrTiO3. Nat. Mater. 2014;13:1085–1090. doi: 10.1038/nmat4107. [PubMed] [Cross Ref]
4. Szot K, Speier W, Bihlmayer G, Waser R. Switching the electrical resistance of individual dislocations in single-crystalline SrTiO3. Nat. Mater. 2006;5:312–320. doi: 10.1038/nmat1614. [PubMed] [Cross Ref]
5. Chen YZ, et al. A high-mobility two-dimensional electron gas at the spinel/perovskite interface of γ-Al2O3/SrTiO3. Nat. Commun. 2013;4:1371. doi: 10.1038/ncomms2394. [PubMed] [Cross Ref]
6. Adler SB. Factors governing oxygen reduction in solid oxide fuel cell cathodes. Chem. Rev. 2004;104:4791–4844. doi: 10.1021/cr020724o. [PubMed] [Cross Ref]
7. Sawa A. Resistive switching in transition metal oxides. Mater. Today. 2008;11:28–36. doi: 10.1016/S1369-7021(08)70119-6. [Cross Ref]
8. Muller DA, Nakagawa N, Ohtomo A, Grazul JL, Hwang HY. Atomic-scale imaging of nanoengineered oxygen vacancy profiles in SrTiO3. Nature. 2004;430:657–661. doi: 10.1038/nature02756. [PubMed] [Cross Ref]
9. Maier J. Nanoionics: ion transport and electrochemical storage in confined systems. Nat. Mater. 2005;4:805–815. doi: 10.1038/nmat1513. [PubMed] [Cross Ref]
10. Kumar A, et al. Probing surface and bulk electrochemical processes on the LaAlO3-SrTiO3 interface. ACS Nano. 2012;6:3841–3852. doi: 10.1021/nn204960c. [PubMed] [Cross Ref]
11. Bark CW, et al. Switchable induced polarization in LaAlO3/SrTiO3 heterostructures. Nano Lett. 2012;12:1765–1771. doi: 10.1021/nl3001088. [PubMed] [Cross Ref]
12. Bühlmann S, Colla E, Muralt P. Polarization reversal due to charge injection in ferroelectric films. Phys. Rev. B. 2005;72:214120. doi: 10.1103/PhysRevB.72.214120. [Cross Ref]
13. Xie Y, Bell C, Hikita Y, Hwang HY. Tuning the electron gas at an oxide heterointerface via free surface charges. Adv. Mater. 2011;23:1744–1747. doi: 10.1002/adma.201004673. [PubMed] [Cross Ref]
14. Bi F, et al. ‘Water-cycle’ mechanism for writing and erasing nanostructures at the LaAlO3/SrTiO3 interface. Appl. Phys. Lett. 2010;97:173110. doi: 10.1063/1.3506509. [Cross Ref]
15. Kim Y, et al. Mechanical control of electroresistive switching. Nano Lett. 2013;13:4068–4074. doi: 10.1021/nl401411r. [PubMed] [Cross Ref]
16. Sharma P, et al. Mechanical tuning of LaAlO3/SrTiO3 interface conductivity. Nano Lett. 2015;15:3547–3551. doi: 10.1021/acs.nanolett.5b01021. [PubMed] [Cross Ref]
17. Morozovska AN, et al. Thermodynamics of electromechanically coupled mixed ionic-electronic conductors: Deformation potential, Vegard strains, and flexoelectric effect. Phys. Rev. B. 2011;83:195313. doi: 10.1103/PhysRevB.83.195313. [Cross Ref]
18. Kogan SM. Piezoelectric effect during inhomogeneous deformation and acoustic scattering of carriers in crystals. Sov. Phys. Solid State. 1964;5:2069–2070.
19. Yudin PV, Tagantsev AK. Fundamentals of flexoelectricity in solids. Nanotechnology. 2013;24:432001. doi: 10.1088/0957-4484/24/43/432001. [PubMed] [Cross Ref]
20. Catalan G, et al. Flexoelectric rotation of polarization in ferroelectric thin films. Nat. Mater. 2011;10:963–967. doi: 10.1038/nmat3141. [PubMed] [Cross Ref]
21. Lee D, et al. Giant flexoelectric effect in ferroelectric epitaxial thin films. Phys. Rev. Lett. 2011;107:057602. doi: 10.1103/PhysRevLett.107.057602. [PubMed] [Cross Ref]
22. Lu H, et al. Mechanical writing of ferroelectric polarization. Science. 2012;336:59–61. doi: 10.1126/science.1218693. [PubMed] [Cross Ref]
23. Lee D, Yang SM, Yoon J, Noh TW. Flexoelectric rectification of charge transport in strain-graded dielectrics. Nano Lett. 2012;12:6436–6440. doi: 10.1021/nl3038129. [PubMed] [Cross Ref]
24. Zubko P, Catalan G, Buckley A, Welche PRL, Scott JF. Strain-gradient-induced polarization in SrTiO3 single crystal. Phys. Rev. Lett. 2007;99:167601. doi: 10.1103/PhysRevLett.99.167601. [PubMed] [Cross Ref]
25. Stengel M. Surface control of flexoelectricity. Phys. Rev. B. 2014;90:201112. doi: 10.1103/PhysRevB.90.201112. [Cross Ref]
26. Sadewasser, S. & Glatzel, T. (eds). Kelvin Probe Force Microscopy: Measuring and Compensating Electrostatic Forces (Springer, 2012).
27. Andrä M, et al. The influence of the local oxygen vacancy concentration on the piezoresponse of strontium titanate thin films. Nanoscale. 2015;7:14351–14357. doi: 10.1039/C5NR03643G. [PubMed] [Cross Ref]
28. Liscio A, Palermo V, Samorì P. Nanoscale quantitative measurement of the potential of charged nanostructures by electrostatic and kelvin probe force microscopy: unraveling electronic processes in complex materials. Acc. Chem. Res. 2010;43:541–550. doi: 10.1021/ar900247p. [PubMed] [Cross Ref]
29. Nielsen DA, Popok VN, Pedersen K. Modelling and experimental verification of tip-induced polarization in Kelvin probe force microscopy measurements on dielectric surfaces. J. Appl. Phys. 2015;118:195301. doi: 10.1063/1.4935811. [Cross Ref]
30. Kim Y, et al. Ionically-mediated electromechanical hysteresis in transition metal oxides. ACS Nano. 2012;6:7026–7033. doi: 10.1021/nn3020757. [PubMed] [Cross Ref]
31. Balke N, et al. Exploring local electrostatic effects with scanning probe microscopy: implications for piezoresponse force microscopy and triboelectricity. ACS Nano. 2014;8:10229–10236. doi: 10.1021/nn505176a. [PubMed] [Cross Ref]
32. Bieger T, Maier J, Waser R. Kinetics of oxygen incorporation in SrTiO3 (Fe-doped): an optical investigation. Sens. Actuators B Chem. 1992;7:763–768. doi: 10.1016/0925-4005(92)80400-R. [Cross Ref]
33. Schie M, Marchewka A, Müller T, De Souza Ra, Waser R. Molecular dynamics simulations of oxygen vacancy diffusion in SrTiO3. J. Phys. Condens. Matter. 2012;24:485002. doi: 10.1088/0953-8984/24/48/485002. [PubMed] [Cross Ref]
34. Paladino AE. Oxidation kinetics of single-crystal SrTiO3. J. Am. Ceram. Soc. 1965;48:476–478. doi: 10.1111/j.1151-2916.1965.tb14803.x. [Cross Ref]
35. Waser R. Bulk conductivity and defect chemistry of acceptor-doped strontium titanate in the quenched state. J. Am. Ceram. Soc. 1991;74:1934–1940. doi: 10.1111/j.1151-2916.1991.tb07812.x. [Cross Ref]
36. Shin C, Yoo H. Al-doped SrTiO3: Part II, unusual thermodynamic factor and chemical diffusivity. Solid State Ion. 2007;178:1089–1094. doi: 10.1016/j.ssi.2007.05.009. [Cross Ref]
37. Kumar A, Ciucci F, Morozovska AN, Kalinin SV, Jesse S. Measuring oxygen reduction/evolution reactions on the nanoscale. Nat. Chem. 2011;3:707–713. doi: 10.1038/nchem.1112. [PubMed] [Cross Ref]
38. Morozovska AN, et al. Flexocoupling impact on size effects of piezoresponse and conductance in mixed-type ferroelectric semiconductors under applied pressure. Phys. Rev. B. 2016;94:174101. doi: 10.1103/PhysRevB.94.174101. [Cross Ref]
39. Freedman DA, Roundy D, Arias TA. Elastic effects of vacancies in strontium titanate: Short- and long-range strain fields, elastic dipole tensors, and chemical strain. Phys. Rev. B. 2009;80:064108. doi: 10.1103/PhysRevB.80.064108. [Cross Ref]
40. Gu Y, Hong Z, Britson J, Chen L-Q. Nanoscale mechanical switching of ferroelectric polarization via flexoelectricity. Appl. Phys. Lett. 2015;106:022904. doi: 10.1063/1.4905837. [Cross Ref]
41. Cao Y, Shen J, Randall C, Chen L-Q. Effect of ferroelectric polarization on ionic transport and resistance degradation in BaTiO3 by phase-field approach. J. Am. Ceram. Soc. 2014;97:3568–3575. doi: 10.1111/jace.13162. [Cross Ref]
42. Wang J-J, et al. Defect chemistry and resistance degradation in Fe-doped SrTiO3 single crystal. Acta Mater. 2016;108:229–240. doi: 10.1016/j.actamat.2016.02.022. [Cross Ref]
43. Cheng G, et al. Electron pairing without superconductivity. Nature. 2015;521:196–199. doi: 10.1038/nature14398. [PubMed] [Cross Ref]
44. Tagantsev AK. Piezoelectricity and flexoelectricity in crystalline dielectrics. Phys. Rev. B. 1986;34:5883–5889. doi: 10.1103/PhysRevB.34.5883. [PubMed] [Cross Ref]
45. Zubko P, Catalan G, Tagantsev AK. Flexoelectric effect in solids. Annu. Rev. Mater. Res. 2013;43:387–421. doi: 10.1146/annurev-matsci-071312-121634. [Cross Ref]
46. Arrhenius SA. Über die dissociationswärme und den einflusß der temperatur auf den dissociationsgrad der elektrolyte. Z. Physik. Chem. 1889;4:96–116.
47. Bednyakov PS, Sluka T, Tagantsev AK, Damjanovic D, Setter N. Formation of charged ferroelectric domain walls with controlled periodicity. Sci. Rep. 2015;5:15819. doi: 10.1038/srep15819. [PMC free article] [PubMed] [Cross Ref]
48. Oh YS, Luo X, Huang F-T, Wang Y, Cheong S-W. Experimental demonstration of hybrid improper ferroelectricity and the presence of abundant charged walls in (Ca,Sr)3Ti2O7 crystals. Nat. Mater. 2015;14:407–413. doi: 10.1038/nmat4168. [PubMed] [Cross Ref]
49. Seidel J, et al. Electronic properties of isosymmetric phase boundaries in highly strained Ca-doped BiFeO3. Adv. Mater. 2014;26:4376–4380. doi: 10.1002/adma.201400557. [PubMed] [Cross Ref]
50. Cao Y, Shen J, Randall C, Chen L-Q. Effect of multi-domain structure on ionic transport, electrostatics, and current evolution in BaTiO3 ferroelectric capacitor. Acta Mater. 2016;112:224–230. doi: 10.1016/j.actamat.2016.04.027. [Cross Ref]
51. Chen L-Q, Shen J. Applications of semi-implicit Fourier-spectral method to phase field equations. Comput. Phys. Commun. 1998;108:147–158. doi: 10.1016/S0010-4655(97)00115-X. [Cross Ref]
52. Li YL, et al. Phase transitions and domain structures in strained pseudocubic (100) SrTiO3 thin films. Phys. Rev. B. 2006;73:184112. doi: 10.1103/PhysRevB.73.184112. [Cross Ref]
53. Moos R, Hardtl KH. Defect chemistry of donor-doped and undoped strontium titanate ceramics between 1000° and 1400°C. J. Am. Ceram. Soc. 1997;80:2549–2562. doi: 10.1111/j.1151-2916.1997.tb03157.x. [Cross Ref]
54. Li YL, Hu SY, Liu ZK, Chen L-Q. Effect of substrate constraint on the stability and evolution of ferroelectric domain structures in thin films. Acta Mater. 2002;50:395–411. doi: 10.1016/S1359-6454(01)00360-3. [Cross Ref]
55. Li YL, Hu SY, Liu ZK, Chen L-Q. Effect of electrical boundary conditions on ferroelectric domain structures in thin films. Appl. Phys. Lett. 2002;81:427–429. doi: 10.1063/1.1492025. [Cross Ref]

Articles from Nature Communications are provided here courtesy of Nature Publishing Group