Phase-field simulations of vortex chirality manipulation in ferroelectric thin films

Chirality, referring to left- and right-handed phenomena, enables diverse applications in optics, sensing, and stereochemistry. Recently, nanoscale textures in ferroelectrics, including polar vortices and skyrmions, have been shown to exhibit chirality and are promising for robust, high-density information storage. BiFeO₃ (BFO), a prototypical rhombohedral ferroelectric with eight <111> polarization variants and 180°, 109°, and 71° domain walls (DWs), serves as an ideal platform to study chiral textures. Prior work has suggested and observed chiral DWs arising from depolarization fields and strain inhomogeneity. Theory indicates curled electric fields can control vortex chirality via toroidal moment switching, while experiments demonstrate creation/annihilation of vortex pairs by local radial electric fields, engineered geometries, and mechanical fields. Flexoelectricity significantly influences polarization switching, DW formation, and vortex morphology, and can be quantified via phase-field methods. Despite progress, deterministic control and reversible switching of ferroelectric vortex chirality in thin films remain challenging. This study uses phase-field simulations to elucidate formation mechanisms and demonstrate controllable manipulation of chiral vortices in BFO thin films under local surface charge or electric field, focusing on how initial bi-domain arrangements and external field direction set vorticity and polarity, respectively, and on the stability and reversibility of chirality.

• Chirality in ferroelectric textures has been observed, including chiral DWs and polar skyrmions/vortices; DW chirality linked to depolarization fields and strain inhomogeneity.
• Curled electric fields can switch toroidal moments to control vortex chirality in nanoparticles; polarization–strain coupling promotes vortex formation.
• In thin films, local radial electric fields and compositional/shape engineering enable creation and annihilation of vortex pairs. Mechanical fields and interface strain can also control vortex morphology and promote polarization rotation to form opposite vortex pairs.
• Flexoelectricity strongly affects polarization switching and DW structures; flexoelectric coupling has been quantified via phase-field modeling and ML, and flexo-sensitive polarization vortices have been reported. Vortex evolution is sensitive to flexoelectricity, which impacts morphology but its role in setting chirality needs clarification.
• Deterministic chirality control in ferroelectrics has been explored in designed nanostructures and magnetic analogs, but robust, reversible chirality switching in ferroelectric thin films with retained topology after stimulus removal remains an open goal that this work addresses.

• Framework: Time-dependent Ginzburg–Landau (TDGL) phase-field model describing polarization evolution Pi(r,t) along x=[100], y=[010], z=[001]: ∂Pi/∂t = −L δftotal/δPi.
• Free energy density: ftotal = fLandau + fgradient + felectric + felastic + fflexo. • Landau energy: sixth-order polynomial in P1,P2,P3 with coefficients α1, α11, α12, α111, α112, α123; α1 = (T−T0)/(2ε0C0). • Gradient energy: Gijkl Pi,j Pk,l simplified with G11, G12, G44 terms accounting for symmetric and antisymmetric polarization gradients. • Flexoelectric energy: fflexo = fijk Pi,j εkl with cubic-symmetry coefficients f11, f12, f44; driving field Eflex derived from ∂fflexo/∂ε producing strain-gradient-coupled terms. • Elastic energy: felastic = Cijkl εij εkl with eigenstrain εij0 = Qijkl Pk Pl (electrostriction). Three elastic constants C11, C12, C44 for cubic symmetry; electrostrictive Q11, Q12, Q44. • Electrostatic energy: felectric = −(1/2) ε0 P·E − (1/2) κi E·E with isotropic background dielectric κ11=κ22=κ33=45; E = −∇φ, φ from Poisson equation ∇²φ = −ρi/(ε0 κi).
• Boundary conditions: • Electrostatic: Short-circuit to avoid depolarization and bound charge screening effects. Local surface charge density σsurf applied on the top surface as free charge boundary: ρiz = nf = σsurf at the surface; ρ = 0 inside film. Surface potential distributions also used to impose local radial electric fields (square-shaped region) for comparison with experiments. • Mechanical: Top surface stress-free; bottom constrained by substrate with equal biaxial misfit strain εxx=εyy=0. Mechanical equilibrium ∂σij/∂xj=0; top surface satisfies ∂σi3/∂x3 = nf = 0.
• Geometries and discretization: • Film thickness: 10 nm (n = 10Δx), grid spacing Δx = 1 nm. • Simulation cells: single-domain 64×64×16, bi-domain 128×128×16, tri-domain 240×240×16 (in units of Δx).
• External stimuli: • Local square-shaped positive/negative surface charge σsurf applied; also local square-shaped surface potential to generate radial electric fields. Reported σsurf values include normalized units (σsurf* from 0 to 10) with conversion σactual = 0.16 C/m² × σsurf*; example images at σsurf*=0.5,1.5,2.5,5,7.5,10 correspond to 0.08, 0.24, 0.40, 0.80, 1.20, 1.60 C/m². • Electric field E* = Ex/0.192 MV/cm used in alternating-field protocols. • Charge region width W varied (40, 60, 80 nm) to assess stability.
• Material parameters: BFO coefficients provided (Table 1): Landau αi, αij, αijk; gradient G11=0.6, G12=0, G44=0.3 (C−2 m² N); elastic c11=2.280×10¹¹, c12=1.250×10¹¹, c44=6.500×10¹⁰ N/m²; electrostrictive Q11=0.032, Q12=−0.016, Q44=0.02; flexoelectric f11=2.5 V, f22=2.5 V, f44=0.05 V; background dielectric κ=45; spontaneous polarization P0=0.52 C/m²; temperature T=300 K.
• Analysis: Domain structures, energy components (Landau, gradient, electrostatic, elastic, surface), vortex core trajectories, and elastic strain fields (ε11, ε22, ε33, ε23, ε13, ε12) evaluated; shear strain maps (ε23, ε13) correlated with vorticity; effects of flexoelectric coupling assessed.

• Single-domain response to local surface charge: • Vortices can be generated from an initial single domain by applying local positive/negative surface charge; chirality is random due to degenerate energy states. • Under positive σsurf (downward built-in field), downward Pz switching in the charged region forms, with upward Pz domains around to minimize stray fields; straight four-quadrant topology evolves to a curved-wall vortex. • Switching pathways: 71° switching dominates over 109° and 180°, yielding a higher percentage of R domains from 71° transitions (e.g., R1+→R3−), consistent with lower energy barriers. • Charge-density dependence: As σsurf increases, domain ratios change. Example: at σ=0.08 C/m² (σ*=0.5), R3− dominates; increasing σ grows R1−. At σ≈0.80 C/m² (σ*=5), a vortex with roughly balanced four R variants forms. Vortex core moves toward the film center as σ increases. • Energy evolution: Increasing σsurf decreases electrostatic energy due to built-in field contributions; more switching increases DW density, raising Landau energy. Curved DWs minimize free and electrostatic energies, stabilizing the vortex in a metastable configuration.
• Bi-domain control of vorticity and chirality: • Domain walls break radial symmetry and deterministically set vorticity. Initial bi-domain orientation determines vortex chirality; reversing the arrangement reverses chirality. • For DWs along y-axis: R3+|R3− → LH (CCW); R3−|R3+ → RH (CW); R3+|R2− → LH; R2−|R3+ → LH. • For DWs along x-axis: R3+/R3− → RH; R3−/R3+ → LH; R3+/R2− → RH; R2−/R3+ → LH. • Equivalent behavior obtained using local square-shaped surface potential (radial electric fields) instead of surface charge.
• Strain–vorticity correlation: • Near vortex cores, ε11 and ε22 are tensile, ε33 is compressive. Shear strains ε23 and ε13 exhibit opposite extrema on either side of the core. The relative directions/signs of ε23 and ε13 correlate with vorticity: opposite directions correspond to CCW or CW as observed; identical directions on both sides correspond to CCW in the studied cases.
• Flexoelectric effect: • Including flexoelectric coupling slightly increases strain gradients but does not change the vortex chirality for a given initial bi-domain (e.g., R3+|R3− remains LH).
• Reversible chirality switching and stability: • Alternating external stimuli (surface charge polarity or electric field sign) reversibly switch vortex chirality with high reproducibility by flipping out-of-plane polarity while preserving in-plane vorticity. • After removing the external field, the topological vortex is retained; morphology shrinks slightly but chirality remains, aided by 71° orientation relationships between inner vortex and surrounding outer domains. • Stability depends on charge-region size: vortices formed with larger W (e.g., 80 nm) remain stable post-removal, whereas smaller regions (W = 40, 60 nm) yield unstable vortices that merge into striped domains.

The simulations show that deterministic control of ferroelectric vortex chirality in BFO thin films can be achieved by leveraging initial bi-domain configurations and external fields. The domain wall orientation and the relative polarization directions across the DW set the in-plane vorticity and thus the handedness of the vortex; reversing the bi-domain arrangement reverses chirality. External field polarity controls the out-of-plane polarization, enabling reversible switching of chirality when combined with a fixed vorticity. The coupling between polarization and elastic strain explains the observed correlations: tensile in-plane and compressive out-of-plane strains near the core accompany continuous polarization rotation, while shear strain patterns (ε23, ε13) encode the vorticity sign. Flexoelectricity influences morphology via enhanced strain gradients but does not alter chirality selection given an initial domain architecture. The retention of the topological vortex after field removal, contingent on sufficient charged-region size, highlights nonvolatile behavior and potential fatigue tolerance critical for memory devices. Charge density governs switching pathways and domain ratios, with 71° switching favored, shaping the resulting vortex morphology and core position. These insights bridge domain engineering (via DWs) with field-driven control to realize robust, switchable chiral textures in ferroelectric films.

This work establishes a deterministic and reversible route to manipulate vortex chirality in BiFeO₃ thin films via phase-field simulations. Key contributions include: (1) elucidating how initial bi-domain arrangements set in-plane vorticity and thus chirality; (2) demonstrating that external field polarity switches out-of-plane polarization to reversibly toggle chirality while maintaining vorticity; (3) revealing the dominant role of 71° switching and the dependence of vortex formation and core position on surface charge density; (4) correlating shear strain signatures with vorticity; and (5) showing post-stimulus retention of topological vortices, with stability governed by the charged-region size. These findings provide theoretical guidance for engineering nonvolatile, fatigue-tolerant chiral vortex states for low-power, high-density ferroelectric memory. Future directions include experimental validation of chirality control protocols in BFO and other ferroelectrics, optimization of device geometries and charge/field delivery to enhance stability at smaller scales, exploration of dynamic switching speeds and endurance, and extension to multilayer or heterostructure systems where interfacial coupling may further tailor vorticity and chirality.

• The study is computational; experimental factors such as defects, leakage, and realistic electrode/ionic environments may influence switching thresholds and stability.
• Stability after field removal depends strongly on the size of the charged region; small regions (W ≤ 60 nm) lead to vortex merging, which may constrain device scaling.
• Material parameters (e.g., flexoelectric coefficients) are based on estimates; variations could affect morphology and thresholds.
• Short-circuit electrostatic boundary conditions and idealized mechanical constraints may not capture all experimental depolarization and substrate effects.