# Phase-Field Study of Nanocavity-Assisted Mechanical Switching in PbTiO 3 Thin Films

Kévin Alhada-Lahbabi, Damien Deleruyelle, Brice Gautier

► To cite this version:

Kévin Alhada-Lahbabi, Damien Deleruyelle, Brice Gautier. Phase-Field Study of Nanocavity-Assisted Mechanical Switching in PbTiO 3 Thin Films. Advanced Electronic Materials, 2024, 10 (5), 10.1002/aelm.202300744. hal-04814170

HAL Id: hal-04814170

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

Submitted on 2 May 2025

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

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

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

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

# Phase-Field Study of Nanocavity-Assisted Mechanical Switching in PbTiO₃ Thin Films

Kévin Alhada-Lahbabi,* Damien Deleruyelle, and Brice Gautier

Ferroelectric thin films hold significant promise in various nanoelectronic applications, demanding precise control over their domain structures. While electrical field-driven polarization switching is currently employed, it often leads to undesirable side effects. In contrast, mechanical switching offers a voltage-free alternative but faces challenges in thicker films. Recent breakthroughs have demonstrated stable mechanical switching in films up to 200 nm thick, attributed to the presence of nanocavities. These nanoscale voids are believed to facilitate domain transitions, serving as essential pinning centers. In this study, mechanical domain switching in thick ferroelectric films is investigated using phase-field modeling, with a specific focus on evaluating the influence of nanocavities on domain stability. The effects of cavity parameters (size, depth, and dielectric properties) on mechanical switching stability under various applied pressures are systematically examined. The findings reveal the intricate interplay between these factors and outline the conditions for stable mechanical switching. Furthermore, phase-field simulations are employed to showcase the energetic mechanisms governing nanocavity-assisted mechanical switching, while also highlighting the pivotal role of these defects as pinning centers. This investigation elucidates the nanocavity-assisted mechanical control of polarization and the potential for optimizing thin film design through nanocavity engineering, thus enabling mechanical switching across substantial film thicknesses.

## 1. Introduction

Ferroelectric materials have demonstrated remarkable potential in an extensive range of nanoelectronic device applications,[1] ranging from nonvolatile memories[2,3] to neuromorphic computing.[4,5] These applications heavily rely on the precise manipulation of ferroelectric domain structures, making a

comprehensive understanding of ferroelectric switching mechanisms an imperative prerequisite. Currently, electrical control of polarization has become a well-established phenomenon in ferroelectric thin films.[2,6,7] Nonetheless, undesirable side effects may arise from using electric fields such as charge injection and breakdown[8,9] or ionic migration,[10] which can hinder the efficacy of domain engineering. To overcome these drawbacks, alternative switching approaches have gained significant attention.[11–13]

Mechanically induced switching has emerged as a promising, voltage-free strategy for locally controlling polarization in ferroelectric materials.[14–24] This approach gained recognition following the pioneering work by Lu et al., who achieved the first polarization reversal in a 4.8-nm-thick BaTiO₃ (BTO) thin film using mechanical means.[17] Subsequently, significant research has been conducted to deterministically control domains through pressure-induced switching, employing atomic force microscopy (AFM) tips to apply force across a diverse range of ferroelectric materials,

including BiFeO₃ (BFO),[14,23] Pb(Zrₓ, Ti₁₋ₓ)O₃ (PZT),[16,18,21] and the recently discovered HfₓZr₁₋ₓO₂ (HZO).[19]

In the realm of ferroelectric thin films, the reversal of polarization is typically ascribed to the flexoelectric field, which arises due to the substantial stress gradient generated beneath the nanoindenter.[17,19,25,26] This effect results in a downward-oriented electric field, enabling a 180-degree unidirectional domain switching that transitions from upward-oriented to downward-oriented polarization. While this effect is typically associated with unidirectional switching, it's noteworthy that recent experimental and theoretical works, have indicated the potential for bidirectional switching,[23,24] thereby potentially expanding the applications of mechanical switching-based devices. However, it's essential to note that the flexoelectric field's influence is limited to a very shallow film thickness, typically around 10 nm.[17,19,25,26] Unfortunately, this limitation hinders the achievement of ferroelectric switching in thicker films, which is crucial for numerous potential commercial applications. The necessity for thicker films arises from the current fabrication processes, aiming to ensure high reliability and performance.

K. Alhada-Lahbabi, D. Deleruyelle, B. Gautier
Univ Lyon, INSA Lyon, CNRS, Ecole Centrale de Lyon, Université Claude Bernard Lyon 1
CPE Lyon, INL, UMR5270, 69621 Villeurbanne, France
E-mail: kevin.alhada-lahbabi@insa-lyon.fr

The ORCID identification number(s) for the author(s) of this article can be found under https://doi.org/10.1002/aelm.202300744

© 2024 The Authors. Advanced Electronic Materials published by Wiley-VCH GmbH. This is an open access article under the terms of the Creative Commons Attribution License, which permits use, distribution and reproduction in any medium, provided the original work is properly cited.

DOI: 10.1002/aelm.202300744
Despite concerted efforts to achieve mechanical polarization control in thicker films, such occurrences have remained relatively rare in the literature, with only a limited number of reported instances.[14,15,18,21,26] In 2016, the control of electrical and mechanical switching was reported in 70-nm-thick BFO films.[14] In a more recent development, Park et al. showcased the ability to selectively control different types of mechanical switching in a 100-nm-thick BFO sample.[26] By utilizing a trailing flexoelectric field induced by a mechanically loaded scanning probe microscope tip, they successfully opted for either a 71-degree or a 180-degree ferroelectric switching scenario. In the field of thick PZT ferroelectric films, Yuhan et al. demonstrated the mechanical reversal of ferroelectric polarization in a 50-nm-thick layer.[21] This achievement was facilitated by the introduction of ferroelastic a domains with in-plane polarization. The incorporation of these a domains within the matrix constituted of out-of-plane c domains effectively lowered the energy barrier for the 180-degree polarization transition. To the best of our knowledge, the most substantial film thickness successfully subjected to mechanical switching was reported by Gonzalez Casal et al. in PbZr<sub>0.2</sub>Ti<sub>0.8</sub>O<sub>3</sub> films, ranging in thickness from 33 to 200 nm.[18] The ability to manipulate domains within films of such considerable thickness was attributed to the presence of nanocavities observed within the sample bulk, positioned in close proximity to the film surface. These void-like nanoscale defects were believed to function as pinning centers for the domains.[18]

Phase-field simulations have become a commonly employed approach for acquiring a deeper understanding of the mechanics dictating polarization dynamics in ferroelectric materials.[27-30] Multiple phase-field endeavors have been undertaken to clarify the fundamental phenomena related to mechanical switching processes.[16,20,23,25,31-40] Regarding mechanical domain control in thick films, these studies have substantially contributed to advancing our understanding of the energetic mechanisms that promote the stability of polarization reversal across significant film thicknesses.[15,18,21,26] Furthermore, several phase-field analyses have showcased the significant effect of the presence of defects, such as voids and cracks, on the behavior of ferroelectric materials.[41-45]

Here, we present a phase-field investigation dedicated to comprehending the pivotal role played by nanocavities in enabling mechanical domain switching in thicker ferroelectric films. We systematically analyze the impact of various cavity parameters on the stability of mechanical switching, elucidating the complex interplay between applied pressure, cavity size, depth from the film surface, and dielectric properties. These investigations are reported in stability diagrams, providing a clear understanding of the essential conditions for successful nanocavity-assisted mechanical switching. Furthermore, our work aims to comprehend the physical phenomena induced by those nanocavities. As such, the evolution of polarization under the influence of tip-induced pressure is illustrated throughout the simulations, shedding light on the underlying energetic mechanisms governing ferroelectric switching stability. This investigation not only enhances our fundamental understanding but also paves the way for the precise control of polarization within void-assisted switching scenarios. We envision these insights holding promise for the design and optimization of films tailored for mechanical switching applications,

positioning nanocavity-assisted domain switching as a key player in the advancement of ferroelectric thin film technology.

## 2. Phase-Field Method

The phase-field simulations are conducted by extending the method presented in a previous study[46] to include mechanical switching modeling. The microstructure evolution of the ferroelectric polarization during domain formations in two dimensions is described by the time-dependent Landau Ginzburg (TDGL) equation[27,28]:

$$ \frac{\partial P_i(\mathbf{r}, t)}{\partial t} = -L \frac{\delta \psi}{\delta P_i(\mathbf{r}, t)}, \quad (i = 1, 2) \qquad (1) $$

In this equation, $P_i(\mathbf{r}, t)$ signifies the spontaneous polarization with two components $P = (P_1, P_2) = (P_x, P_y)$, $L$ is a kinetic coefficient, $\psi$ is the total free energy and $\mathbf{r} = (x, y)$ denotes the spatial vector in two dimensions. Here, $\mathbf{x}$ and $\mathbf{y}$ denote the vertical and horizontal axes, respectively. The total free energy $\psi$ encompasses contributions from bulk ($\psi_{\text{bulk}}$), gradient ($\psi_{\text{grad}}$), electric ($\psi_{\text{elec}}$), elastic ($\psi_{\text{elas}}$), and flexoelectric ($\psi_{\text{flexo}}$) free energy density:

$$ \psi = \iint_V \left( \psi_{\text{bulk}}(P_i) + \psi_{\text{grad}}(P_{i,j}) + \psi_{\text{elec}}(P_i, E_i) + \psi_{\text{elas}}(P_i, \epsilon_{ij}) \right. \\ \left. \qquad + \psi_{\text{flexo}}(P_i, P_{i,j}, \epsilon_{ij}, \epsilon_{ij,k}) \right) \delta V \qquad (2) $$

where $P_i, E_i$, and $\epsilon_{ij}$ represent the polarization, electric field, and strain components, respectively. $P_{i,j}$ and $\epsilon_{ij,k}$ denote the polarization and strain gradient, with the notation $P_{i,j}$ defined as $P_{i,j} = \frac{\partial P_i}{\partial x_j}$. In two dimensions, the bulk energy is described by Ref. [27]:

$$ \psi_{\text{bulk}}(\mathbf{P}) = \alpha_1(P_1^2 + P_2^2) + \alpha_{11}(P_1^4 + P_2^4) + \alpha_{111}(P_1^6 + P_2^6) + \alpha_{12}P_1^2 P_2^2 \\ + \alpha_{112}(P_1^4 P_2^2 + P_1^2 P_2^4) \qquad (3) $$

where $\alpha_{11}, \alpha_{12}, \alpha_{111}$, and $\alpha_{112}$ are the Landau coefficients, and $\alpha_1 = \frac{T-T_0}{2\epsilon C}$ refers to the dielectric permittivity $\epsilon$, the Curie temperature $T_0$ and the Curie constant $C$. The energy contribution of the domain walls is described by the gradient energy, which in a cubic system is calculated by Ref. [27]:

$$ \psi_{\text{grad}}(P_{i,j}) = \frac{1}{2} G_{11} (P_{1,1}^2 + P_{2,2}^2) + G_{12} P_{1,1} P_{2,2} + \frac{1}{2} G_{44} (P_{2,1} + P_{1,2})^2 \\ + \frac{1}{2} G'_{44} (P_{2,1} - P_{1,2})^2 \qquad (4) $$

where $G_{11}, G_{12}$, and $G_{44}$ are the gradient energy coefficients. The electrostatic energy is described by the following equation:

$$ \psi_{\text{elec}}(P_i, E_i) = -P_i E_i - \frac{1}{2} \epsilon_0 \epsilon_r E_i E_j \qquad (5) $$

here, $E_i$ represents the electric field in the thin film, $\epsilon_0$ is the vacuum dielectric permittivity, and $\epsilon_r$ is the relative dielectric permittivity. The electric field $E_i = -\frac{\partial V}{\partial x_i}$ is derived from the
electrostatic potential V, which is obtained by solving the electrostatic equilibrium[47]:

$$ \frac{\partial}{\partial x_i} [\epsilon_0 \epsilon_r(\mathbf{r}) \frac{\partial V(\mathbf{r})}{\partial x_j}] = \frac{\partial P_i(\mathbf{r})}{\partial x_i} - \rho(\mathbf{r}) \quad (6) $$

in which $\rho$ represents the electric charge density, and $\frac{\partial P_i}{\partial x_i}$ denotes the polarization charge density. The electrostatic equilibrium is solved using the fast Fourier transform spectral method described in Ref. [47]. The ferroelectric layer is assumed to be under short-circuit boundary conditions, where both the bottom electrode is grounded and the top surface is under complete charge screening[24] ($V_{\text{bot}} = V_{\text{top}} = 0$ V). Additionally, as the thin film is assumed to be infinite in the horizontal direction, periodic boundary conditions are applied along with the in-plane direction.

The elastic energy density is given by Ref. [27, 30]:

$$ \psi_{\text{elas}}(P_i, \epsilon_{ij}) = \frac{1}{2} C_{ijkl} (\epsilon_{ij}(\mathbf{r}) - \epsilon_{ij}^0(\mathbf{r})) (\epsilon_{kl}(\mathbf{r}) - \epsilon_{kl}^0(\mathbf{r})) \quad (7) $$

where $C_{ijkl}$ represents the elastic stiffness tensor, $\epsilon_{ij}(\mathbf{r})$ represents the total strain, and $\epsilon_{ij}^0(\mathbf{r}) = Q_{ijkl} P_k(\mathbf{r}) P_l(\mathbf{r})$ is the electrostrictive strain caused by the polarization, where $Q_{ijkl}$ is the electrostrictive tensor. The total strain $\epsilon_{ij}(\mathbf{r})$ comprises both the homogeneous strain $\epsilon_{ij}^s$ and the heterogeneous strain $\epsilon_{ij}(\mathbf{r}) = \delta\epsilon_{ij}(\mathbf{r}) = \frac{1}{2}(\frac{\partial u_i(\mathbf{r})}{\partial x_i} + \frac{\partial u_j(\mathbf{r})}{\partial x_j})$, where $u_i$ represents the components of the mechanical displacement field.

In the case of ferroelectric materials, the mechanical equilibrium equation $\sigma_{ij,j} = 0$ can be written as Ref. [27, 30]:

$$ C_{ijkl} \frac{\partial^2 u_k(\mathbf{r})}{\partial x_i \partial x_j} = C_{ijkl} \frac{\partial \epsilon_{kl}^0(\mathbf{r})}{\partial x_j} \quad (8) $$

In this study, thin film boundary conditions are employed, as in other related investigations,[25,30–32,34] and the displacements induced by the tip are directly incorporated into the boundary conditions of the mechanical equilibrium.[31,32,34] At the top surface of the thin films of thickness $h$, beneath the tip ($x = h$), the stress is prescribed, and the mechanical boundary condition is given as:

$$ \sigma_{ij} n_i = \tau_i^{\text{tip}} \quad (9) $$

here, $n_i$ represents the normal vector to the top surface, and $\tau_i$ denotes the components of the surface traction induced by the tip load. The contact load is modeled following the Hertz model, assuming an isotropic linear material,[25,32,34,39] where the vertical induced pressure is given by:

$$ \tau_1^{\text{tip}}|_{x=h} = \begin{cases} -\frac{3p}{2\pi a^2} \sqrt{1 - \left(\frac{r}{a}\right)^2} & \text{if } r < a \\ 0 & \text{if } r \ge a \end{cases} \quad (10) $$

where $p$ represents the induced-force-load, $a$ is the contact area radius, and $r$ denotes the distance from the tip center. Additional information about the modeling of the contact between the tip

<table><caption>**Table 1.** Energy coefficients and input parameters of the phase-field simulation.</caption><thead><tr><th>Landau coefficients</th><th>Elastic/Electrostrictive</th><th>Simulation Parameters</th></tr></thead><tbody><tr><td>α<sub>1</sub> (C<sup>-2</sup>.m<sup>2</sup>.N)</td><td>C<sub>11</sub> (GPa)</td><td>P<sub>0</sub> (Cm<sup>-2</sup>)</td><td>0.757</td></tr><tr><td>α<sub>11</sub> (C<sup>-4</sup>.m<sup>6</sup>.N)</td><td>C<sub>12</sub> (GPa)</td><td>ε<sub>r</sub> (a.u.)</td><td>50</td></tr><tr><td>α<sub>12</sub> (C<sup>-4</sup>.m<sup>6</sup>.N)</td><td>C<sub>44</sub> (GPa)</td><td>1111.1</td><td>f<sub>11</sub>, f<sub>12</sub>, f<sub>44</sub> (V)</td><td>10, 0, 0</td></tr><tr><td>α<sub>111</sub> (C<sup>-6</sup>.m<sup>10</sup>.N)</td><td>Q<sub>11</sub> (C<sup>-2</sup>.m<sup>4</sup>)</td><td>0.089</td><td>R<sub>tip</sub> (nm)</td><td>10</td></tr><tr><td>α<sub>112</sub> (C<sup>-6</sup>.m<sup>10</sup>.N)</td><td>Q<sub>12</sub> (C<sup>-2</sup>.m<sup>4</sup>)</td><td>-0.026</td><td>p (μN)</td><td>1–4</td></tr><tr><td>α<sub>123</sub> (C<sup>-6</sup>.m<sup>10</sup>.N)</td><td>Q<sub>44</sub> (C<sup>-2</sup>.m<sup>4</sup>)</td><td>0.03375</td><td>t<sub>FE</sub> (nm)</td><td>5–64</td></tr></tbody></table>

and the sample can be found in the Supporting Information. On the top surface surrounding the tip or in the absence of load, stress-free boundary conditions are applied. Moreover, the mechanical perturbation induced by the tip is assumed to vanish in the substrate. Therefore, a nil displacement is prescribed far from the interface of the film with the substrate.[30–32] Finally, periodic boundary conditions on the mechanical displacements are applied along the in-plane direction. The mechanical equilibrium was solved using the finite element methods (FEM) through the simulations, using the SfePy python library.[48]

In the presence of important strain gradients beneath the tip, the ferroelectric polarization is subject to influences not solely attributable to the piezoelectric effect but also arising from the contribution of the flexoelectric effect:

$$ P_i = d_{ijk} \epsilon_{jk} + \mu_{ijkl} \frac{\partial \epsilon_{kl}}{\partial x_j} \quad (11) $$

where $d_{ijk}$ represents the piezoelectric tensors, and $\mu_{ijkl}$ denotes the flexoelectric tensor. The flexoelectric energy, which is mainly induced by the strain gradient, can be expressed as Ref. [25]:

$$ \psi_{\text{flexo}} = \frac{f_{ijkl}}{2} \left( \frac{\partial P_k}{\partial x_i} \epsilon_{ij} - \frac{\partial \epsilon_{ij}}{\partial x_i} P_k \right) = \frac{F_{ijkl}}{2} \left( \frac{\partial P_k}{\partial x_i} \sigma_{ij} - \frac{\partial \sigma_{ij}}{\partial x_i} P_k \right) \quad (12) $$

where $f_{ijkl}$ and $F_{ijkl}$ are the flexocoupling coefficient tensors, which are related by the relations $f_{ijkl} = C_{ijklm} F_{mnkl}$[39] and $\mu_{ijkl} = \epsilon_0 \epsilon_{ijkl} f_{ijkl}$.[15]

In this paper, a PbTiO₃ thin film is modeled, with material parameters chosen as in Refs. [39, 46, 49] and listed in Table 1. The 2D system is discretized on a grid of size $h \Delta x \times 128 \Delta y$, where $t_{FE} = h \Delta x$ denotes the film thickness and $\Delta x$ and $\Delta y$ represent the grid spacing. The film thickness varies within a range of 5–64 nm across the different cases explored throughout the manuscript. Further details on the simulation parameters are provided in Supporting Information. An illustration of the mechanical displacement and stress distribution in a 20-nm-thick film under an applied load of 2 µN is presented in Figure 1. In this specific example, the influence of ferroelectric polarization on mechanical displacements and stresses within the layer was intentionally disregarded, as the objective was to isolate the deformation induced by the tip load. However, it is crucial to note that the polarization contribution is thoroughly addressed in the remaining sections of the manuscript. Figure 1a,b depicts the substantial displacements beneath the AFM tip. Notably, vertical displacements up to several angstroms[25] are attained just under the tip. Due to the extremely small tip radius of just a few nanometers,
<figure><img src="image_3.png" /><figcaption>**Figure 1.** Mechanical Displacements and Stress Fields in Response to an Induced Tip-Pressure. The contribution of ferroelectric polarization to mechanical displacements and stress distribution is not considered in this figure, prioritizing clarity in depicting the tip-induced perturbation. a) Spatial distribution of mechanical displacements within a 20-nm-thick film subjected to an induced tip-pressure of 2 µN applied over a contact area with a tip radius of *a* = 10 nm. Displacements are magnified tenfold for clarity compared to the colorbar representation. The substrate, where zero displacement is imposed far from the interface, is not depicted. b) Spatial distribution of mechanical displacements (*u*<sub>x</sub>, *u*<sub>y</sub>) within the thin film close to the tip. c) Spatial distribution of stress components (*σ*<sub>11</sub>, *σ*<sub>12</sub>) within the thin film close to the tip.</figcaption></figure>

significant stress values are observed in close proximity to the tip, as shown in Figure 1c, reaching several GPa of intensity.

**3. Results and Discussion**

**3.1. Mechanical Switching in Nanocavity-Free Thin Films**

In this section, we offer an overview of the mechanisms govern-
ing mechanical switching in ferroelectric thin films in the ab-

sence of nanocavities. We begin by examining the response of
a 10-nm-thick ferroelectric film to mechanical switching under
a 4 µN tip load. **Figure 2** depicts the polarization evolution at
various time steps across the simulation. Immediately after tip
loading (Figure 2a), the vertical polarization component is sup-
pressed beneath the tip.[³⁹] The polarization transition evolves
into a 90-degree switching, transforming the initial c⁺ domains
into a ferroelectric domains, as seen in Figure 2b. This switch-
ing is facilitated by the substantial elastic energy generated by

<figure><img src="image_4.png" /><figcaption>**Figure 2.** Dynamical Evolution of Polarization and Stable Domain Switching in a 10-nm-Thick Film Under an Applied Pressure of 4 µN. a) At *t* = 30Δ*t*, immediately following the application of pressure, the vertical component of the polarization is suppressed beneath the tip. b) At *t* = 1000Δ*t*, the polarization undergoes a transition from *c*⁺ to a ferroelectric domains throughout the entire film thickness. c) At *t* = 1050Δ*t*, shortly after the removal of the tip, a 90-degree switching from *a* to *c*⁻ domains occurs. d) At *t* = 1500Δ*t*, the mechanically switched *c*⁻ domain remains stable upon tip removing.</figcaption></figure>
<figure><img src="image_5.png" /><figcaption>**Figure 3.** Dynamical Evolution of Polarization and Unstable Domain Switching in a 20-nm-Thick Film Under an Applied Pressure of 4 µN. a) At *t* = 30Δ*t*, following the pressure application, the vertical component of the polarization is suppressed beneath the tip. b) At *t* = 1000Δ*t*, the polarization undergoes a 90-degree switching. Unlike the 10-nm-thick film, this transition does not extend to the bottom electrode. c) At *t* = 1050Δ*t*, the polarization switches back to a vertical *c*<sup>+</sup> domain. d) At *t* = 1500Δ*t*, the polarization has reverted to its initial *c*<sup>+</sup> domain structure.</figcaption></figure>

compressive stress beneath the tip, favoring in-plane polariza-
tion. Upon tip removal (Figure 2c), the polarization undergoes a
second 90-degree switching, transitioning from **a** to **c**<sup>-</sup> domains.
This shift is primarily driven by the persistent influence of misfit
strain, compelling the polarization to adopt a vertical orientation.
Remarkably, for the 10-nm-thick film, the initial **c**<sup>+</sup> → **a** 90-degree
switching extends throughout the entire film thickness. Impor-
tantly, in this configuration, the polarization does not induce a
depolarizing field when the **c**<sup>-</sup> domains arise, thanks to the bot-
tom electrode's role as a charge compensation region. Overall,
this process results in a 180-degree switch from the initial state,
aided by the absence of depolarizing charge due to the clamped
domain on the bottom electrode. The mechanically switched **c**<sup>-</sup>
domain then remains stable long after tip removal (Figure 2d).
These dynamics demonstrate the feasibility of stable mechanical
switching in thin layers, consistent with previous literature.[39]

In contrast, previous reports[39] have indicated that such a
mechanism is unlikely to occur in thicker films, as exempli-
fied by the 20-nm-thick film in **Figure 3**. Here, the initial **c**<sup>+</sup> →
**a** 90-degree switching does not penetrate the full film depth,
remaining confined to approximately 10 nm beneath the tip
(Figure 3a,b). Upon tip removal, the polarization reverts to a ver-
tical **c**<sup>+</sup> domain, driven by the strong upward-oriented depolar-
izing field (Figure 3c,d). Indeed, the nucleation of **c**<sup>-</sup> domains
would generate significant depolarizing charges, hindering me-
chanical switching. These observations underscore the critical
role of film thickness in determining the stability of mechanical
domain switching.

In order to provide a comprehensive overview of domain
switching behavior across varying film thicknesses and applied
pressures, we present a representation of the criteria for me-
chanical switching stability in **Figure 4**. This diagram illustrates
the stability of mechanically switched domains within thin films,
spanning thicknesses from 5 to 25 nm, and subject to tip loads
varying between 1 and 4 µN. Our findings align with established
trends observed in related studies,[39] consistently indicating that
mechanical switching tends to become unachievable in thicker
films. Furthermore, these results highlight the importance of

an adequate load, as mechanical pressure can significantly af-
fect switching stability for certain film thicknesses. For exam-
ple, the 5-nm-thick film becomes unstable to mechanical switch-
ing when the tip-induced force exceeds 3 µN. In contrast, a 15-
nm-thick film requires a force above 1 µN in order to initiate a
90-degree switching mechanism that propagates throughout the
thickness and allows stable mechanical switching. Further illus-
trations of those switching scenarios are provided in Figures S1–
S4 (Supporting Information). Specifically, the physical mecha-
nisms underlying the observed instability in the 5-nm-thick film
under a 4 µN mechanical load are thoroughly explained in the

<figure><img src="image_6.png" /><figcaption>**Figure 4.** Switchable and Nonswitchable Domains in Thin Films without Nanocavities. A dot symbolizes stable c⁻ domain switching that persists after the removal of the tip. In contrast, a cross depicts unstable domain switching, characterized by the polarization reverting to an upward-oriented c⁺ domain after the tip's removal. Boundaries in the stability map are interpolated between points, with simulations solely conducted at cross and dot points.</figcaption></figure>
<figure><img src="image_7.png" /><figcaption>**Figure 5.** Dynamical Evolution of Polarization and Stable Domain Switching within an 8-nm-Radius Nanocavity Under an Applied Pressure of 2 µN. a) At $t = 30\Delta t$, the polarization suppression beneath the tip reaches the void. b) At $t = 1000\Delta t$, the 90-degree polarization switching extends to the top of the nanocavity. c) At $t = 1050\Delta t$, the **a** domain undergoes a 90-degree switching clamping onto the cavity. d) At $t = 1500\Delta t$, the **c**<sup>-</sup> domain is stabilized, anchored on the top surface of the void.</figcaption></figure>

Supporting Information. Slight deviations from the threshold values observed in the literature[39] are noticeable, although the general trends remain consistent. These discrepancies can be attributed to variations in material parameters, mechanical boundary conditions, and the expression of the mechanical stress induced by the AFM tip.

Notably, the stability map boundaries in Figure 4 are interpo-
lated from simulations conducted exclusively at cross and dot
points. A finer mesh in the studied parameters space would
yield smoother boundaries, improving the accuracy of phase di-
agrams. However, it is worth stressing that achieving this refine-
ment entails substantial computational costs, given the complex-
ity and resource demands inherent in phase-field simulations.

**3.2. Nanocavity-Assisted Mechanical Switching**

3.2.1. Influence of Cavity Size and Depth on Mechanical Switching Stability

Nanocavities stand as promising factors that could potentially
contribute to the stabilization of mechanically switched domains
within thicker ferroelectric films.[18] In this section, we delve into
the compelling potential of nanocavities, when suitably incorpo-
rated, to enhance and stabilize mechanical domain switching in
thick films well beyond 20 nm in thickness. Through a com-
prehensive analysis, we aim to shed light on the intricate in-
terplay between nanocavities and mechanical domain switching,
offering insights into the mechanisms that drive their stabiliz-
ing effects. We focused on spherical nanocavities to stay consis-
tent with experimental observations.[18] In the phase-field simu-
lations, nanocavities were treated as dielectric voids where the fer-
roelectric polarization was assumed to be nil (P = 0). For simplic-

ity, their dielectric and elastic properties were initially assumed
identical to the ones of the PbTiO₃ thin film in this section. Later,
cavities with lower permittivity will be considered. The polariza-
tion behavior within the cavity, along with the electrostatic poten-
tial, mechanical displacements, and stress distributions, is fur-
ther elucidated in Figure S5 (Supporting Information).

To illustrate the interplay between nanocavities and mechan-
ical domain switching, we provide in **Figure 5a** detailed repre-
sentation of the evolution dynamics of polarization in a 64-nm-
thick film subjected to an applied load of 2 µN. A nanocavity,
featuring an 8-nm radius and centered 32 nm beneath the sur-
face, is introduced into the film. Following the application of tip
load, the polarization initially undergoes the same $c^+ \rightarrow a$ tran-
sition as described in the anterior sections (Figure 5a). The 90-
degree polarization switching process extends over an approx-
imate length of 20 nm, encompassing the upper region of the
nanocavity (Figure 5b). Upon removal of the tip, the polarization
must undergo a second 90-degree transition to attain a vertical
state, driven by the strong compression misfit. However, in the
presence of the cavity, the depolarizing field experienced if under-
going a $a \rightarrow c^-$ transition by the switched domain is significantly
reduced. Importantly, this clamping effect counters the head-to-
head configuration that often impedes mechanical switching in
thick films without nanocavities, thereby creating energetically
favorable conditions for the nucleation of a fully switched $c^-$
domain (Figure 5c). Finally, we observe the stabilization of the
mechanically switched $c^-$ domain, which remains pinned onto
the nanocavity (Figure 5d). However, in contrast to the defect-
free films, this mechanically switched domain remains confined
to the space between the top of the nanocavity and the film's
surface.

In the following section, we examine the influence of
nanocavity size, location, and applied pressure on stability. Our
<figure><img src="image_8.png" /><figcaption>**Figure 6.** Switchable and Nonswitchable Domains in Thin Films with Nanocavities. Diagrams showcasing the stability of mechanically switched domains within 64-nm-thick thin films under tip pressures ranging from 1 to 4 µN, highlighting the influence of nanocavity size and location on stability (a), as the influence of cavity dielectric permittivity (b). A dot signifies stable c⁻ domain switching that endures even after the removal of the tip, whereas a cross represents unstable domain switching, characterized by the polarization reverting to an upward-oriented c⁺ domain following the tip's removal. Boundaries in the phase diagrams are interpolated between points, with simulations solely conducted at cross and dot points.</figcaption></figure>

investigation focuses on 64-nm-thick thin films subjected to mechanical loads ranging from 1 to 4 µN. To understand the role of nanocavities, we vary the nanocavity radius from 4 to 8 nm, positioning their centers at distances from the film surface spanning 10 to 30 nm. The results are summarized in **Figure 6**, which provides a comprehensive overview of the stability of mechanically switched domains within these thin films. One of our key findings underscores the pivotal role of nanocavity size in determining the stability of mechanically switched domains. A sufficiently large nanocavity is a prerequisite for domain clamping, establishing energetically favorable conditions for domain nucleation and stability. This favorable outcome is achieved by reducing the depolarizing field, facilitated by the cavity's role as a pinning center. This effect is more pronounced in the case of larger cavities, which provide a larger clamping area. Additionally, for the reversed domain to remain stable on the cavity, it must be of sufficient size to avoid negative impacts from gradient energy. This condition is more easily met with larger cavities that have a greater radius. Conversely, smaller nanocavities fail to provide the necessary energetic conditions for domain clamping, rendering it more challenging for the switched domain to maintain stability. This trend is evidenced in Figure 6a, where the largest 8-nm-radius cavity exhibits stability in mechanical switching across various depths and pressures, while mechanical switching is consis-

tently unachievable in the smallest 4-nm-radius cavity. Between these two size extremes, the switching effectiveness becomes contingent on the applied pressure and the depth at which the nanocavity is situated.

In the context of cavity depth, our findings emphasize that the deeper the nanocavity is positioned within the film, the more challenging it becomes to achieve mechanical switching, regardless of other parameters. Similarly to the scenario observed in cavity-free thick films, this phenomenon is attributed to the stress concentration primarily confined around the tip's vicinity. If the domain switching process fails to extend to the location of the nanocavity, it remains unclamped, leading to the generation of strong depolarizing charges and causing the domain to revert to its initial state.

Our observations also reveal that higher applied pressures do not necessarily equate to enhanced switchability. Instead, our simulations unveil a complex and intricate relationship between these three critical parameters. As an illustrative example, the stability of a 5-nm-radius void positioned 10 nm from the film surface is elucidated in Figure 6a. It demonstrates the void's robustness under a 1 µN tip load while highlighting its loss of stability with an increased tip pressure of 2 µN. The polarization-driven switching process is depicted for both cases in **Figure 7**. A noteworthy distinction between the stable (Figure 7a–d) and unstable
<figure><img src="image_9.png" /><figcaption>**Figure 7.** Impact of Mechanical Pressure on Nanocavity-Assisted Domain Switching. a–d) Sequential snapshots illustrating the evolution of the polarization in a stable domain switching scenario for a 5-nm-radius cavity located 10 nm from the surface under a 1 µN tip-load, at a) *t* = 30Δ*t*, b) *t* = 1000Δ*t*, c) *t* = 1050Δ*t*, and d) *t* = 1500Δ*t*. e–h) Sequential snapshots captured at identical times, depicting the evolution of polarization in an unstable domain switching scenario with the same cavity under a 2 µN tip-pressure.</figcaption></figure>

(Figure 7e–h) scenarios arises from the stress distribution, which notably spreads under the cavity in the 2 µN pressure regime. This complexity highlights the need for fine-tuning the applied pressure to optimize stress distribution, elastic/flexoelectric energies, and, ultimately, achieve the desired stability for a particular combination of cavity size and depth.

### 3.2.2. Influence of Cavity Dielectric Permittivity on Mechanical Switching Stability

In the following section of our investigation, we delve into nanocavity dielectric permittivity and its impact on mechanical

domain switching. This exploration is crucial due to the inherent uncertainty surrounding the precise dielectric properties of these void-like defects. While our previous analyses assumed, for simplicity, that nanocavity permittivity matched that of the surrounding PbTiO₃ layer, more realistic physical considerations would suggest it may be significantly lower, possibly aligning more closely with the void's dielectric properties. In Figure 6b, we present the outcomes of our investigation into nanocavities with varying relative permittivity ε_r values, ranging from 1 to 50. These films were subjected to the same tip pressures as previously described, while the nanocavities themselves had radius ranging from 4 to 8 nm and were consistently positioned at a depth of 30 nm from the film's top surface. It's important to note that
<figure><img src="image_10.png" /><figcaption>(c)</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>(e)</figcaption></figure>

<figure><img src="image_12.png" /><figcaption>(g)</figcaption></figure>

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

<figure><img src="image_14.png" /><figcaption>(d)</figcaption></figure>

<figure><img src="image_15.png" /><figcaption>(f)</figcaption></figure>

<figure><img src="image_16.png" /><figcaption>(h)</figcaption></figure>

<figure><img src="image_17.png" /><figcaption>**Figure 8.** Impact of Dielectric Permittivity on Domain Switching under a 1 µN tip-load. a–d) Sequential snapshots illustrating the evolution of the polarization in an unstable domain switching scenario, at a) $t = 30\Delta t$, b) $t = 1000\Delta t$, c) $t = 1050\Delta t$, and d) $t = 1500\Delta t$. The relative permittivity of the 6-nm-radius nanocavity is set 50, aligning with the background permittivity of PbTiO₃. e–h) Sequential snapshots captured at identical times, depicting the evolution of polarization in a stable domain switching scenario facilitated by the same cavity but featuring a reduced relative permittivity of five.</figcaption></figure>

throughout our simulations, we maintained a consistent background relative permittivity of 50 for the ferroelectric PZT layer.

Overall, nanocavities characterized by lower relative permittivity values demonstrate an enhanced propensity for mechanical switchability. This trend remains consistent across various tip loads and cavity sizes, as visually depicted in Figure 6b. Remarkably, in several cases for constant load and cavity size, the emergence of a discernible permittivity threshold becomes apparent. This phenomenon is notably pronounced for the 6-nm-radius cavity, highlighting the delicate interplay between cavity properties required for stable switching at this cavity size.

To gain deeper insights into the underlying physical mechanisms driving domain switching in the presence of low-permittivity voids, we illustrate in Figure 8a comparative analysis of unstable and stable switching scenarios employing a consistent 6-nm-radius nanocavity subjected to a 1 µN indentor load. The key distinction between these cases lies in the relative permittivity of the nanocavity material: it is 50 for the unstable switching (Figure 8a–d) and 5 for the stable switching (Figure 8e–h). As depicted, the critical observation is that lower nanocavity permittivity results in a more pronounced depolarizing field at the cavity boundary. Consequently, polarization must align more rigorously with the cavity’s boundary, aiming to
<figure><img src="image_18.png" /><figcaption>**Figure 9.** Influence of Multiple Cavities on Domain Switching Stability. a–d) Sequential snapshots illustrating the evolution of polarization in a stable domain switching scenario under a 4 µN tip-load at a) *t* = 30Δ*t*, b) *t* = 1000Δ*t*, c) *t* = 1050Δ*t*, and d) *t* = 1500Δ*t*. The interplay resulting from the presence of two 8-nm-radius cavities, located slightly to the right and left of the tip axis, promotes mechanical switching stability upon tip removal.</figcaption></figure>

counteract the polarization flux at the cavity surface. This effect interplays with the compressive epitaxial strain, which inherently favors maintaining vertical polarization, leading to more pronounced horizontalization of the initially upward-oriented c⁺ domain above the nanocavity in the low-permittivity scenario (Figure 8e) than in the high-permittivity case (Figure 8a). This, in turn, promotes the extension of the 90-degree switch induced by the applied pressure, as the presence of strong depolarizing charges created by the original c⁺ domain is significantly reduced. Upon tip removal, this favorable electrostatic landscape further aids the nucleation and growth of the c⁻ domain (Figure 8e–h), a process that would have been hindered in the case of a higher cavity permittivity (Figure 8a–d). To provide a comprehensive overview of the energetic mechanisms occurring during mechanical switching, we elucidate, in the supporting figure Figure S6 (Supporting Information), the evolution of the free energy curves for both stable (ε<sub>r</sub> = 5) and unstable (ε<sub>r</sub> = 50) scenarios.

### 3.2.3. Exploring Multicavity Effects on the Stability of Mechanical Switching

In experimental samples, the presence of multiple nanocavities is commonly observed,¹⁸¹ emphasizing the need to understand the interaction between those distinct cavities and its impact on mechanical switching stability. A thorough analysis of the influence of multiple cavities, accounting for variations in size, location, and permittivity, necessitates extensive and demanding high-throughput computations, and will be the topic of future research efforts.

To offer an initial glimpse into the potential implications of considering multiple cavities, we present a simplified case investigating two closely located cavities, as depicted in Figure 9. In this configuration, two 8-nm-radius cavities (with permittivity ε<sub>r</sub> = 50) are introduced into the film, positioned 20 nm away from the film surface, under a tip force of 4 µN. Importantly, the cav-

ities were intentionally shifted approximately 10 nm to the left and right of the tip axis to explore the interplay of the switched domain between the voids. When the cavities are shifted in this manner from the tip axis, an isolated cavity fails to promote mechanical switching, as illustrated in Figures S7 and S8 (Supporting Information), underscoring the significance of the cavity’s relative location to the tip axis.

Remarkably, the presence of two cavities results in a stable switching state after tip removal, as illustrated in Figure 9. During loading (Figure 9a,b), the polarization undergoes a c⁺ → a transition, leading to the formation of two a domains positioned above the cavities’ top surfaces. Upon tip removal (Figure 9c), each a domain undergoes a transition toward a c⁻ fully switched domain, with each domain stabilized on one of the defects. After relaxation (Figure 9d), two stably mechanically switched domains can be observed, each one clamped on a void surface. This illustrative example, wherein the presence of multiple cavities enhances the stability of mechanical switching, highlights the intricate interplay resulting from diverse distributions of cavities and voids within the film.

## 4. Conclusion

In summary, a phase-field investigation was conducted to elucidate the intricate physical mechanisms governing nanocavity-assisted mechanical switching in ferroelectric thick films. We found that, under specific conditions, nanocavities energetically facilitate the stabilization of switched domains, serving as crucial pinning centers that enhance domain clamping. Our comprehensive findings emphasize the complex interplay between nanocavity size, depth from the film surface, dielectric relative permittivity, and applied tip force. Through systematic phase-field simulations, we have provided a detailed account of how each of these parameters influences mechanical switching stability.

Furthermore, this study also provides a comprehensive overview of the essential conditions required for achieving
nanocavity-assisted mechanical switching in thick films, encapsulated in the form of stability diagrams.

In this study, we conducted simulations on 2D systems, aligning with related research on mechanical switching in thin films.[31,32,35] Given the axial symmetry of the indenter, mechanical field distributions are anticipated to exhibit strong agreement with 3D computations.[32] While this symmetry argument holds for spherical cavities aligned with the tip axis, it becomes invalid for arbitrary cavity shapes and distributions. Extending the model to three dimensions thus becomes necessary for obtaining more quantitative results, and will be addressed in future research endeavors.

Our investigation also sheds light on the underlying energetic effects that govern the stability of ferroelectric switching in the context of nanocavities. Our findings underscored the significant role of the depolarizing energy in hindering stable ferroelectric switching. Importantly, nanocavities have been demonstrated as an effective means to overcome this depolarizing field, enabling efficient switching without being constrained by film thickness.

Nanocavity-engineered ferroelectric films may present exciting potential applications, such as integration into ferroelectric memories. In these films, mechanical switching could enable unidirectional voltage-free operations, while opposite direction switching could be achieved using electrical stresses. Additionally, these films hold promise for nanometric-resolution pressure sensors, utilized, for instance, as impact detectors. The incorporation of nanocavity-assisted switching enables nanoscale impacts to trigger mechanical switching, surpassing a threshold force and generating distinguishable switched domains. Notably, the key advantage of nanocavity-switching over traditional mechanical switching lies in its capability to facilitate mechanical switching in thick films, reaching up to hundreds of nanometers. This characteristic would greatly streamline the integration of mechanical-switching devices into modern fabrication processes. However, it is crucial to acknowledge that practical applications of mechanical switching remain constrained by physical limitations, such as inducing reliable bidirectional switching in a voltage-free manner. Although nanocavity-assisted switching alone may not overcome all these limitations, envisioning future studies that integrate cavities and recent insights in mechanical bidirectional switching[23,24] holds promise for addressing these challenges and unlocking new possibilities in the field.

By enhancing our understanding of the prerequisites for achieving stable void-assisted pressure-induced switching in thick films, this work contributes to the efficient design of thick ferroelectric films that incorporate nanocavities. This advancement paves the way for applications reliant on mechanical switching in thick films, opening new horizons for the utilization of ferroelectric materials in various technological domains.

## Supporting Information

Supporting Information is available from the Wiley Online Library or from the author.

## Conflict of Interest

The authors declare no conflict of interest.

## Data Availability Statement

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

## Keywords

cavities, ferroelectrics, mechanical switching, PbTiO₃, phase-field modeling

Received: October 25, 2023

Revised: January 9, 2024

Published online: January 31, 2024

[1] J. F. Scott, C. A. P. de Araujo, *Science* **1989**, *246*, 1400.

[2] A. Crassous, S. Sluka, A. K. Tagantsev, N. Setter, *Nat. Nanotechnol.* **2015**, *10*, 614.

[3] P. Sharma, D. Sando, Q. Zhang, X. Cheng, S. Prosandeev, R. Bulanadi, S. Prokhorenko, L. Bellaiche, L.-Q. Chen, V. Nagarajan, J. Seidel, *Adv. Funct. Mater.* **2019**, *29*, 1807523.

[4] T. Mikolajick, S. Slesazeck, H. Mulaosmanovic, M. H. Park, S. Fichtner, P. D. Lomenzo, M. Hoffmann, U. Schroeder, *J. Appl. Phys.* **2021**, *129*, 100901.

[5] M. H. Park, Y. H. Lee, H. J. Kim, Y. J. Kim, T. Moon, K. D. Kim, J. Müller, A. Kersch, U. Schroeder, T. Mikolajick, C. S. Hwang, *Adv. Mater.* **2015**, *27*, 1811.

[6] J. Ma, J. Ma, Q. Zhang, R. Peng, J. Wang, C. Liu, M. Wang, N. Li, M. Chen, X. Cheng, P. Gao, L. Gu, L.-Q. Chen, P. Yu, J. Zhang, C.-W. Nan, *Nat. Nanotechnol.* **2018**, *13*, 947.

[7] L. McGilly, P. Yudin, L. Feigl, A. Tagantsev, N. Setter, *Nat. Nanotechnol.* **2015**, *10*, 145.

[8] A. Q. Jiang, Y. Y. Lin, T. A. Tang, *J. Appl. Phys.* **2007**, *102*, 074109.

[9] T. Kawae, Y. Terauchi, H. Tsuda, M. Kumeda, A. Morimoto, *Appl. Phys. Lett.* **2009**, *94*, 112904.

[10] W. Ming, B. Huang, S. Zheng, Y. Bai, J. Wang, J. Wang, J. Li, *Sci. Adv.* **2022**, *8*, eabq1232.

[11] L. Ponet, S. Artyukhin, T. Kain, J. Wettstein, A. Pimenov, A. Shuvaev, X. Wang, S.-W. Cheong, M. Mostovoy, A. Pimenov, *Nature 2022*, *607*, 81.

[12] Y.-D. Liou, Y.-Y. Chiu, R. T. Hart, C.-Y. Kuo, Y.-L. Huang, Y.-C. Wu, R. V. Chopdekar, H.-J. Liu, A. Tanaka, C.-T. Chen, C. F. Chang, L. H. Tjeng, Y. Cao, V. Nagarajan, Y.-H. Chu, Y. C. Chen, J. C. Yang, *Nat. Mater.* **2019**, *18*, 580.

[13] Y. Tian, L. Wei, Q. Zhang, H. Huang, Y. Zhang, H. Zhou, F. Ma, L. Gu, S. Meng, L.-Q. Chen, C.-W. Nan, J. Zhang, *Nat. Commun.* **2018**, *9*, 3809.

[14] L. Chen, Z. Cheng, W. Xu, X. Meng, G. Yuan, J. Liu, *Sci. Rep.* **2016**, *6*, 19092.

[15] Y. Lun, X. Wang, J. Kang, Q. Ren, T. Wang, W. Han, Z. Gao, H. Huang, Y. Chen, L.-Q. Chen, D. Fang, J. Hong, *Adv. Mater.* **2023**, *35*, 2302320.

[16] G. Vats, Ravikant, P. Schoenherr, A. Kumar, J. Seidel, *Adv. Electron. Mater.* **2020**, *6*, 2000523.

[17] H. Lu, C.-W. Bark, D. E. de los Ojos, J. Alcala, C. B. Eom, G. Catalan, A. Gruverman, *Science* **2012**, *336*, 59.

[18] S. Gonzalez Casal, X. Bai, K. Alhada-Lahbabi, B. Canut, B. Vilquin, P. Rojo Romeo, S. Brottet, D. Albertini, D. Deleruyelle, M. Bugnet, I. Canero Infante, B. Gautier, *Adv. Electron. Mater.* **2022**, *8*, 2200077.

[19] Z. Guan, Y.-K. Li, Y.-F. Zhao, Y. Peng, G. Han, N. Zhong, P.-H. Xiang, J.-H. Chu, C.-G. Duan, *Nano Lett.* **2022**, *22*, 4792.

[20] Z. Wen, X. Qiu, C. Li, C. Zheng, X. Ge, A. Li, D. Wu, *Appl. Phys. Lett.* **2014**, *104*, 042907.
[21] G. Yuan, H. Huang, C. Li, D. Liu, Z. Cheng, D. Wu, *Adv. Electron. Mater.* **2020**, *6*, 2000300.

[22] H. Lu, B. Wang, T. Li, A. Lipatov, H. Lee, A. Rajapitamahuni, R. Xu, X. Hong, S. Farokhipoor, L. W. Martin, C.-B. Eom, L.-Q. Chen, A. Sinitskii, A. Gruverman, *Nano Lett.* **2016**, *16*, 6460.

[23] Y. Wang, C. Guo, M. Chen, Y. Liang, H. Han, H. Chen, Y. Lin, D. Yi, H. Huang, C.-W. Nan, J. Ma, *Adv. Funct. Mater.* **2023**, *33*, 2213787.

[24] J. Liu, W. Chen, M. Wu, F. Sun, X. Huang, Y. Zheng, *npj Comput. Mater.* **2022**, *8*, 137.

[25] Y. Gu, Z. Hong, J. Britson, L.-Q. Chen, *Appl. Phys. Lett.* **2015**, *106*, 022904.

[26] S. M. Park, B. Wang, S. Das, S. C. Chae, J. S. Chung, J. G. Yoon, L. Q. Chen, S. M. Yang, T. W. Noh, *Nat. Nanotechnol.* **2018**, *13*, 366.

[27] L.-Q. Chen, *Annu. Rev. Mater. Res.* **2002**, *32*, 113.

[28] J.-J. Wang, B. Wang, L.-Q. Chen, *Annu. Rev. Mater. Res.* **2019**, *49*, 127.

[29] Y. Zhao, *npj Comput. Mater.* **2023**, *9*, 94.

[30] Y. L. Li, S. Y. Hu, Z. K. Liu, L. Q. Chen, *Acta Mater.* **2002**, *50*, 395.

[31] J. Li, W. Xiong, X. Huang, W. Chen, Y. Zheng, *J. Appl. Phys.* **2021**, *129*, 244105.

[32] C. Qi, Y. Jiang, X. Wang, C. S. Lynch, *J. Appl. Phys.* **2022**, *131*, 244101.

[33] L. Jiang, J. Tang, Y. Zhou, Q. Yang, Y. Zhang, L. Guo, X. Zhong, *Comput. Mater. Sci.* **2015**, *108*, 309.

[34] W. Xiong, J. Liu, L. Ma, W. Chen, Y. Zheng, *J. Appl. Phys.* **2020**, *128*, 014102.

[35] C. Qi, X. Wang, *J. Appl. Phys.* **2021**, *129*, 094101.

[36] L. L. Ma, W. J. Chen, B. Wang, W. M. Xiong, Y. Zheng, *J. Phys.: Condens. Matter* **2019**, *32*, 035402.

[37] C. Wang, X. Ke, J. Wang, R. Liang, Z. Luo, Y. Tian, D. Yi, Q. Zhang, J. Wang, X.-f. Han, G. V. Tendeloo, L.-q. Chen, C.-w. Nan, *Nat. Commun.* **2016**, *7*, 10636.

[38] C. Guo, S. Dong, J. Wang, X. Wang, H. Huang, *Adv. Electron. Mater.* **2023**, *9*, 2300193.

[39] Y. Cao, A. Morozovska, S. V. Kalinin, *Phys. Rev. B* **2017**, *96*, 184109.

[40] Y. Cao, S. V. Kalinin, *Phys. Rev. B* **2016**, *94*, 235444.

[41] Y. C. Song, A. K. Soh, Y. Ni, *J. Phys. D: Appl. Phys.* **2007**, *40*, 1175.

[42] A. Abdollahi, I. Arias, *J. Mech. Phys. Solids* **2012**, *60*, 2100.

[43] S. Yuan, W. Chen, L. Ma, Y. Ji, W. Xiong, J. Liu, Y. Liu, B. Wang, Y. Zheng, *Acta Mater.* **2018**, *148*, 330.

[44] P. Fedeli, M. Kamlah, A. Frangi, *Smart Mater. Struct.* **2019**, *28*, 035021.

[45] R. Indergand, D. M. Kochmann, *Appl. Phys. Lett.* **2021**, *119*, 222901.

[46] K. Alhada-Lahbabi, D. Deleruyelle, B. Gautier, *ACS Appl. Electron. Mater.* **2023**, *5*, 3894.

[47] J. J. Wang, X. Q. Ma, Q. Li, J. Britson, L. Q. Chen, *Acta Mater.* **2013**, *61*, 7591.

[48] R. Cimrman, V. Lukeš, E. Rohan, *Adv. Comput. Math.* **2019**, *45*, 1897.

[49] J. Wang, S. Q. Shi, L. Q. Chen, Y. Li, T. Y. Zhang, *Acta Mater.* **2004**, *52*, 749.