=begin = Molecular Dynamics Simulations of Hysteresis Loops for BaTiO3 Ferroelectric Thin-Film Capacitors Using the feram Code Takeshi Nishimatsu http://dx.doi.org/10.1103/PhysRevB.78.104104 http://loto.sourceforge.net/feram/ CLICK anywhere to go to the next page. # I entitled my talk "Molecular Dynamics Simulations of Hysteresis Loops for BaTiO3 Ferroelectric Thin-Film Capacitors Using the feram Code". # This study is done by my original MD code, feram. # You can download the code from here. == How to use Slidy This presentation is written with Slidy http://www.w3.org/Talks/Tools/Slidy/ . Use Firefox web browser. * Advance to next slide with mouse click or [Space] * Move forward/backward between slides with [Left], [Right], [Pg Up] and [Pg Dn] keys * [Home] key for first slide, [End] key for last slide * The [C] key for an automatically generated table of contents (or click on "contents" on the toolbar) * Function [F11] toggles between the full screen mode and the normal mode * [F] key toggles the display of the tool bar below. * [A] key toggles the presentation mode and the printing mode * [S] and [B] keys can control the size of fonts == What is ferroelectric? Fig. 1 figures/ferroelectrics.jpg Properties of ferroelectrics, such as BaTiO$_3$, PbTiO$_3$, KNbO$_3$, LiNbO$_3$, LiTaO$_3$, SrTiO$_3$, Pb(Zr$_{(1-x)}$Ti$_x$)O$_3$ (PZT), etc.. /Fig. # First, let me remind you what ferroelectric is. # As you know, # ferroelectrics such as BaTiO3, PbTiO3, SrTiO3, PZT, etc. # has characteristic properties, for example, # spontaneous polarization # phase transition # piezoelectric effect and # inverse piezoelectric effect. # # You can see both piezoelectric and inverse piezoelectric effects with this demonstrator. # This white part is PZT. # If you apply AC voltage, it can be a speaker. # And it can be a sensor, too. # # Here, I'd like to emphasize that # without long-range dipole-dipole interaction and finite size effects, # you will not see the domain structure like this. # In other words, the domain structure is formed to avoid # depolarization field induced by surface charges. == Non-volatile ferroelectric memory (FeRAM) Fig. 1 figures/twoFeRAM.jpg Schematic device configurations of two types of FeRAM (a) 1-transistor and 1-capacitor (1T1C) type and (b) 1-transistor (1T) type. (c) and (d) are their circuit diagrams, respectively. /Fig. # Such ferroelectrics are currently used in Non-volatile ferroelectric memory, FeRAM. # They are schematic device configurations and their circuit diagrams. # This one called 1-transistor and 1-capacitor type and # this DRAM like configuration called 1-transistor type. # Both of them has ferroelectric thin-film capacitors. == Ferroelectric thin-film capacitors (background of this study) * Various applications of ferroelectric thin-films: * multilayer capacitors * nonvolatile FeRAMs * nanoactuators * Down-sizing of FeRAMs (nano-capacitors of ferroelectric thin films) is highly demanded. * Epitaxially constraint ferroelectric films made by PLD. * Effect of imperfect screening of electrodes is NOT well understood. * [M. Dawber et al.: J. Phys. Condes. Matter 15, L393 (2003)] * Fatigue causes damages in polarization behaviors of ferroelectric capacitors. Fig. 1 figures/J.F.Scott.Fig.1.1.ae.jpg After [J. F. Scott: Ferroelectric Memories (Springer, 2000)]. /Fig. # I explain the background of this study with following two slides. # There are various applications of ferroelectric thin films. == Dynamics, nanosize effects and temperature dependences * Dynamics of ferroelectric thin-film capacitors: * Hysteresis loop * Polarization switching * Dynamics of domain wall * Their nanosize effects and temperature dependences remain poorly known. * Experimentally, in situ observations are difficult. * Theoretically, * The long-range Coulomb interaction limits the size and time of molecular-dynamics (MD) simulations. * How to include depolarization fields caused by imperfect electrodes was unclear. Fig. 1 figures/Fong.jpg After [Fong et al.: Science 304, 1650 (2004)]. Ultra thin films ($< 4$ layers) may not have spontaneous polarization. /Fig. # There are many dynamical behaviors in ferroelectric thin-film capacitors. # For example, Hysteresis loop, polarization switching, and # dynamics of domain wall. # # For example, this is ... == Purpose of this study * Development of fast molecular dynamics (MD) code which can simulate ferroelectric thin-film capacitors for a realistic system size (up to 100 nm) and a realistic time span ($>$ 1 ns). * Clarify the effect of dead layers between ferroelectrics and electrodes. * Predict the effect of the compressive strain arising from epitaxial constraints. $\to$ We develop "feram" code. Fig. 1 ../logo.jpg /Fig. # Again, one of the purposes of this simulation program is to study # the thermal properties # of ferroelectric thin film and its capacitors. # # Ferroelectric thin-film capacitors are applied, in for example, non-volatile FeRAM. # And their down-sizing are highly demanded. # # Fatigue of electrodes of capacitors is a current big problem # # But their Hysteresis, polarization switching, dynamics of domain wall # and their nano-size effects are not well known. # # This is an example of experimental study of thin films. # Fong et al. studied PbTiO3 thin films' striped ferroelectric domain structure # with strong X-ray. They found thickness dependence of Tc. # # Basically, experimental {\em in situ} observations of such dynamics are difficult. # # But you can do simulation with feram. == Features of feram program * Scientific features * Molecular dynamics (MD) simulation with first-principles-based effective Hamiltonian * $AB$O3-type perovskite ferroelectrics and relaxors * Uses supercell; ($AB$O3)$_N$, $N=32\times 32\times 512$ $\Rightarrow$ 2,621,440 atoms * Coarse-graining; reduction of the number of degree of freedom * Long-range dipole-dipole interaction is treated in reciprocal-space; k-locality of the force matrix * Applications: * Phase transition of bulk ferroelectrics * Capacitor, ferroelectric thin film is sandwiched between short-circuited electrodes * Technical features * Fast * FFTW library version 3 http://www.FFTW.org/ * Parallelized with OpenMP http://www.OpenMP.org/ * Multi-platform * Linux PC * HITACHI SR11000 Super Technical Server * SONY PLAYSTATION3 (ongoing) * Object oriented programming (OOP) in Fortran 2003 * GNU autotools http://www.gnu.org/software/autoconf/ * Free software (GPLv3) http://loto.sf.net/feram/ # Other features of feram program are: # Its Molecular dynamics (MD) simulation uses first-principles-based effective Hamiltonian. # It can simulate ABO3-type perovskite ferroelectrics and relaxors # It can use big size of supercell up to 32x32x512. It corresponds to 2,621,440 atoms. # It can do such a huge calculation because it employs some kind of Coarse-graining that # reduces the number of degree of freedom of the system. # and because Long-range dipole-dipole interaction is treated # in reciprocal-space using k-locality of the force matrix # # Using this program, you can simulate # Phase transition of bulk ferroelectrics and properties of # capacitors in which ferroelectric thin film is sandwiched between short-circuited electrodes. # Today, I will mainly show you our study of ferroelectric thin films. # # This program is fast, because it uses a fast FFT library and # it is parallelized with OpenMP. # I usually use normal Linux personal computers. # But you can also use feram on some supercomputers. # You can freely download feram from this site. == Soft-mode displacements of perovskite type ferroelectrics ABO3 # [King-Smith and Vanderbilt: Phys. Rev. B 49 5828 (1994)] Parametrization of coupling between soft-mode displacement and strains. Fig. 1 figures/perovskite-contours.jpg Coupling between atomic displacements and strains. After [T. Hashimoto, T. Nishimatsu et al.: Jpn. J. Appl. Phys. 43, 6785 (2004)]. /Fig. Eq. 1 u_{\alpha}^2 = \{v_{\alpha}^A\}^2 + \{v_{\alpha}^B\}^2 + \{v_{\alpha}^{\rm O_{I}}\}^2 + \{v_{\alpha}^{\rm O_{II}}\}^2 + \{v_{\alpha}^{\rm O_{III}}\}^2 /Eq. # The origin of electric polarization of perovskite type ferroelectrics ABO3 is its Soft-mode displacement. # And the atomic displacement is strongly coupled with strains. # This is the cubic high temperature crystal structure of ABO3. # Unit cell has 5 atoms. So atomic displacements spans a 3x5=15-dimensional space. # This is a 2-dimensional subspace of the 15-dimensional space of atomic displacements. # In the 2-dimensional subspace, potential surface is plotted. # If there is strain in crystal, the potential minimum goes far from its origin and get deeper. # # With first-principles calculations, you can investigate such potential surface and then construct the effective Hamiltonian. # In this paper, I described how to trace this valley line automatically with first-principles calculations. # This is very useful to construct the effective Hamiltonian. == Potential surface of the zone-center distortion for BaTiO3 Fig. 1 figures/comp-bto-thick.jpg Calculated potential surfaces of the zone-center distortion for BaTiO3 with various $E_{\rm xc}$. /Fig. == Coarse-graining * Real perovskite-like system has $15N+6$ degree of freedom * $N$ unit cells in a supercell * 5 atoms per unit cell * Each atom can move along 3 directions * 6 components of strain * Simplified model has $6N+6$ degree of freedom * Define 1 dipole $Z^*\bm{u}(\bm{R})$ per unit cell * 1 acoustic displacement (Inhomogeneous strain) $\bm{w}(\bm{R})$ per unit cell Fig. 1 figures/CoarseGraining.jpg /Fig. # In the Hamiltonian, Coarse-graining is employed. # In other words, the number of degree of freedom is reduced. # # Real perovskite-like system has 15N+6 degree of freedom, # where N is the number of unit cells in supercell, # because each unit cell has 5 atoms, each atom can move along 3 directions, # and there are 6 components of homogeneous strain. # # For Fast and Simple simulation, we reduce the number of degree of freedom to 6N+6, # by replacing 5 atoms to 2 vectors, 1 dipole Z^*u(R) and 1 acoustic displacement # per unit cell. # == First-principles effective Hamiltonian for a supercell Supercell of ($AB$O3)$_N$ with $\{\bm{u}(\bm{R})\}$ and $\{\bm{w}(\bm{R})\}$: Eq. 1 H^{\rm eff} = \frac{M^*_{\rm dipole}}{2} \sum_{\bm{R},\alpha}\dot{u}_\alpha^2(\bm{R}) + \frac{M^*_{\rm acoustic}}{2}\sum_{\bm{R},\alpha}\dot{w}_\alpha^2(\bm{R}) + V^{\rm self}(\{\bm{u}\})+V^{\rm dpl}(\{\bm{u}\})+V^{\rm short}(\{\bm{u}\}) + V^{\rm elas,\,homo}(\eta_1,...,\eta_6)+V^{\rm elas,\,inho}(\{\bm{w}\}) + V^{\rm coup,\,homo}(\{\bm{u}\}, \eta_1,...,\eta_6)+V^{\rm coup,\,inho}(\{\bm{u}\}, \{\bm{w}\}) - Z^*\sum_{\bm{R}} \mathscr{E} \cdot \bm{u}(\bm{R}) /Eq. Fig. 1 figures/arrows.jpg Snapshot of the supercell. Strengths and directions of $\{\bm{u}(\bm{R})\}$ are indicated with arrows. red and blue colors represent $+z$-polarized and $-z$-polarized sites. /Fig. # Using that parameters, we can construct an effective Hamiltonian for a supercell with # N unit cells of ABO3. # # This is a snapshot of an MD simulation. u(R)s are indicated with arrows. # We do time evolution of this Hamiltonian of this supercell. == Self potential of a local dipole Eq. 1 V^{\rm self}(\bm{u})=\kappa_2 u^2 + \alpha u^4 + \gamma(u_{y}^2u_{z}^2+u_{z}^2u_{x}^2+u_{x}^2u_{y}^2) /Eq. where Eq. 1 u^2=u_{x}^2 + u_{y}^2 + u_{z}^2. /Eq. Fig. 1 figures/fourth.jpg Schematic illustration of self potential of a local dipole. Double well is made by the fourth order terms of u. /Fig. # I shall introduce you some terms of the effective Hamiltonian. # The self potential of a local dipole or On-site local potential has fourth order terms of u. # It makes double well like this for u. == Dipole-dipole interaction In the effective Hamiltonian, dipole-dipole interactions are separated into the long-rage term and the short-range term. Long-range: Eq. 1 V^{\rm dpl}(\{\bm{u}\}) = \frac{1}{2}\sum_{i=1}^N \sum_\alpha \sum_{j=1}^N \sum_\beta u_\alpha(\bm{R}_i) \Phi_{\alpha\beta}(\bm{R}_{ij}) u_\beta(\bm{R}_j), /Eq. Eq. 2 \Phi_{\alpha\beta}(\bm{R}_{ij})= \frac{Z^{*2}}{\epsilon_\infty}\sum_{\bm{n}}{}' \frac{\delta_{\alpha\beta} - 3(\widehat{\bm{R}_{ij}+\bm{n}})_\alpha(\widehat{\bm{R}_{ij}+\bm{n}})_\beta}{|{\bm{R}_{ij}+\bm{n}}|^3}, /Eq. $\bm{n}$ is the supercell lattice vector: Eq. 3 n_\alpha = \cdots, -2{L_\alpha} a_0, -{L_\alpha} a_0, 0, {L_\alpha} a_0, 2{L_\alpha} a_0,\cdots /Eq. Short-range: Eq. 4 V^{\rm shrt}(\{\bm{u}\})= \frac{1}{2}\sum_{i=1}^N \sum_\alpha \sum_j^{\rm 3nn} \sum_\beta u_\alpha(\bm{R}_i) \,J_{ij,\alpha\beta}\, u_\beta(\bm{R}_j) /Eq. $J_{ij,\alpha\beta} = J_k$: Short-range interaction matrix ($k=1,\cdots,7$) # In the effective Hamiltonian, dipole-dipole interactions are separated into the long-rage term and the short-range term. # Both has quadratic form. # But long-range one must be summed-up for very far periodic supercells # Short range one is local. # == Forces on dipoles are calculated in the reciprocal space Fig. 1 figures/flow.jpg Simplified flow chart for calculating forces on dipoles. Fast Fourier transformation (FFT) and inverse FFT (IFFT) enable the rapid calculation of long-range dipole-dipole interactions. You can reduce computational time from $O(N^2)$ to $O(N \log N)$. /Fig. # Time consuming part in simulation is the calculation of forces of dipoles because # dipole-dipole interactions are long-range. ..... You can reduce computational time # from O(N^2) to O(NlogN). == Long-range dipole-dipole interaction Fig. 1 figures/bulk32x32x32.dipole-dipole-long.interaction.jpg Long-range interaction energy per unit cell is plotted in the first Brillouin zone for the simple cubic crystal. Only with dipole-dipole interaction, the minimum is at the M point, i.e., the checker-board-type antiferroelectric order is most stable. /Fig. == Long-range + short-range interaction Fig. 1 figures/bulk32x32x32.dipole-dipole-long+short.interaction.jpg By introducing short-range interaction, the minimum come to the $\Gamma$ point, i.e., the ferroelectric order is most stable. /Fig. == First Brillouin zone for simple cubic crystals Fig. 1 figures/sc.BZ.jpg First Brillouin zone for simple cubic crystals. /Fig. == LO-TO splitting LO-TO splitting is the ion version of the plasma frequency. Fig. 1 figures/oscillation.jpg Schematic illustration of longitudinal oscillation of a dielectric thin film. Polarization $P_z=Z^* u_z / \Omega$ perpendicular to the film causes surface charges $\pm P_z$. The surface charges result in depolarization field $E_{\rm d}=-4\pi P_z$ in the film. This depolarization field causes restoring force $F=-4\pi Pz Z^* u_z=-4\pi Z^{*2} u_z^2 / \Omega$ on $u_z$. /Fig. #The surface charges result in depolarization field $E_{\rm d}=-\frac{4\pi}{\epsilon_\infty} P_z$ in the film. This depolarization field causes restoring force $F=-\frac{4\pi}{\epsilon_\infty} Pz Z^* u_z=-\frac{4\pi}{\epsilon_\infty} Z^{*2} u_z^2 / \Omega$ on $u_z$. == Boundary condition for capacitors Fig. 1 figures/capacitor.jpg Ferroelectrics are sandwiched between short-circuited perfect electrodes (a) and imperfect electrodes (b). In (a) depolarization field $E_{\rm d}=0$, but in (b) $E_{\rm d}=-4\pi P_z \frac{d}{l+d}$. /Fig. # This slide show you how to simulate thin films # sandwiched between short-circuited electrodes. # As we learned when we were under-graduate students, # An electrode act like a mirror for charges. # So two parallel electrodes act like two parallel mirrors. # You have infinite mirror images beyond the mirrors. # Such configuration can simulate with doubled-size supercell # which have TWO l thickness. # # moreover, if you introduce a gap d in the doubled-size supercell, # you can automatically simulate damaged or imperfect electrodes. == Parameters for BaTiO3 Hamiltonian (input file) * [King-Smith and Vanderbilt: Phys. Rev. B 49, 5828 (1994)] * [Zhong, Vanderbilt, and Rabe: Phys. Rev. B 52, 6301 (1995)] #--- Method, Temperature, and mass --------------- method = 'md' # 'md' - Molecular Dynamics with Nose-Poincare thermostat (Canonical ensemble) # 'lf' - Leap Frog (Microcanonical ensemble) # 'hl' - Hysteresis Loop # 'mc' - Monte Carlo (not implemented yet) GPa = -5.0 kelvin = 300 mass_amu = 39.0 # Required for MD Q_Nose = 1.0 #--- System geometry ----------------------------- bulk_or_film = 'bulk' L = 32 32 32 a0 = 3.94 lattice constant a0 [Angstrom] #--- Time step ----------------------------------- dt = 0.002 [ps] n_thermalize = 40000 n_average = 10000 n_coord_freq = 5000 Write a snapshot to the disk every 5000 steps #--- On-site (Polynomial of order 4) ------------- P_kappa2 = 5.502 [eV/Angstrom^2] # P_4(u) = kappa2*u^2 + alpha*u^4 P_alpha = 110.4 [eV/Angstrom^4] # + gamma*(u_y*u_z+u_z*u_x+u_x*u_y), P_gamma = -163.1 [eV/Angstrom^4] # where u^2 = u_x^2 + u_y^2 + u_z^2 #--- Inter-site ---------------------------------- j = -2.648 3.894 0.898 -0.789 0.562 0.358 0.179 j(i) [eV/Angstrom^2] #--- Elastic Constants --------------------------- B11 = 126. B12 = 44.9 B44 = 50.3 [eV] #--- Elastic Coupling ---------------------------- B1xx = -211. [eV/Angstrom^2] B1yy = -19.3 [eV/Angstrom^2] B4yz = -7.75 [eV/Angstrom^2] #--- Dipole -------------------------------------- init_dipo_avg = 0.00 0.00 0.00 [Angstrom] # Average of initial dipole displacements init_dipo_dev = 0.02 0.02 0.02 [Angstrom] # Deviation of initial dipole displacements Z_star = 9.956 epsilon_inf = 5.24 == Time evolution and fluctuation Fig. 1 figures/plot.epit-16x16-15-01-5.0GPa-0.01-cooling.jpg Energy versus time-step plot of a cooling simulation. You can see Nose-Poincare Hamiltonian remains zero. /Fig. # This figure shows the time steps of simulations. # First, the system exchanges energy between the # one particle heat bath, then it reaches the equilibrium. # You can see Nose-Poincare Hamiltonian remains zero. == Results: Monte Carlo (MC) vs. Molecular Dynamics (MD) feram is a molecular dynamics (MD) program for bulk and thin-film ferroelectrics. Fig. 1 figures/MC-vs-MD.jpg Simulations of bulk BaTiO3 with MC and MD. Horizontal axis is temperature. Vertical axis is strain. While Monte Carlo steps do NOT have physical meaning, MD simulation steps have a physical meaning of the time evolution. Cooling and heating MD simulation can show hysteresis in temperature clearly. /Fig. Fig. 1 figures/BaTiO3-abc-T-dependence.jpg Experimentally observed one after [H. E. Kay and P. Vousden: Philos. Mag. 40, 1019 (1949)]. The relatively weak first-order nature of the cubic-to-tetragonal phase transition can be compared with the first-order tetragonal-to-orthorhombic and orthorhombic-to-rhombohedral phase transitions. /Fig. # feram is a molecular dynamics program for ferroelectrics. # # This shows simulations of bulk BaTiO3 with MC and MD. # Horizontal axis is temperature. Vertical axis is strain. # Compared to old Monte Carlo simulations, MD simulation steps have a physical # meaning of the time evolution. So, for example, Cooling and heating MD simulation # can show hysteresis in temperature, # as you may often see in textbooks of ferroelectrics. == Dielectric constants can be calculated from fluctuations Fig. 1 figures/strain-susceptibility.jpg Temperature dependence of three components of strain and relative dielectric constant of bulk BaTiO3. Three phase transitions are clearly seen. Dielectric constants can be calculated from fluctuations of dipoles. /Fig. == Heating-up and cooling-down simulations of ferroelectric capacitors Animations of horizontal slices of heating-up and cooling-down simulations for BaTiO3 thin-film capacitors with short-circuited electrodes under 1% in-plane biaxial compressive strain. The +z-polarized and −z-polarized sites are denoted by red open squares and blue filled squares, respectively. Used supercell sizes are Lx×Ly×Lz = 40×40×2(l+d) . * (a) l=15, d=1 * heating-up epit-40x40-L015-D1-5.0GPa-heating-animation.gif * cooling-down epit-40x40-L015-D1-5.0GPa-cooling-animation.gif * (b) l=31, d=1 * heating-up epit-40x40-L031-D1-5.0GPa-heating-animation.gif * cooling-down epit-40x40-L031-D1-5.0GPa-cooling-animation.gif * (c) l=127, d=1 * heating-up epit-40x40-L127-D1-5.0GPa-heating-animation.gif * cooling-down epit-40x40-L127-D1-5.0GPa-cooling-animation.gif * (d) l=255, d=1 * heating-up (coming soon) * cooling-down epit-40x40-L255-D1-5.0GPa-cooling-animation.gif * (e) l=32, d=0 * heating-up epit-40x40-L032-D0-5.0GPa-heating-animation.gif * cooling-down epit-40x40-L032-D0-5.0GPa-cooling-animation.gif == Single domain structures and striped domain structures Fig. 1 figures/epit-32x32-L015-D1-5.0GPa-100K-three.jpg Snapshots at $T=100$ K in heating-up ((a)) and cooling-down ((b) and (c)) simulations of ferroelectric thin-film capacitors of $l=15$ with $d=1$. (a) and (b) are horizontal slices. (c) is a vertical cross section. /Fig. == Thinner film has finer domain structure to avoid stronger depolarization field Fig. 1 figures/slice32x32.jpg Horizontal slices of domain structures in capacitors of epitaxial BaTiO3 thin-film with imperfect electrodes $32\times 32 \times (l=7, 15, 31, 63, 127 and 255, d=1$). Thinner film has finer domain structure to avoid stronger depolarization field $E_{\rm d}=-4\pi P_z \frac{d}{l+d}$. /Fig. # Here I show thickness dependence of striped domain structure. # We found that thinner film has finer domain # structure to avoid stronger depolarization field. == Temperature dependence of polarizations of thin-film capacitors Fig. 1 figures/epitDEN.jpg /Fig. # These complicated figures are the calculated results for the thin-film # capacitors. Horizontal axes are temperature and vertical axes are # averaged dipole length u-z. # Red solid line can be considered as z-component of # polarization for heating-up simulation. Blue one is that for cooling. # # You can find the thickness dependent transition for heating-up simulations. # but no transition for cooling-down simulations. # # It can be understood there is thermally unhoppable bistability between # uniformly out-of-plane-polarized ferroelectric structure and # striped ferroelectric domain structure. # # In the next slide, I show the snapshot of this point structure. == Uniformly polarized structure and striped domain structure Fig. 1 figures/epit-32x32-L015-D1-5.0GPa-100K.jpg Almost uniformly polarized structure ((a) and (b)) and striped domain structure ((c) and (d)) of thin-film BaTiO3 capacitors of $32\times 32 \times (l=15, d=1)$. There is thermally unhoppable bistability between these two structure. /Fig. # these are the horizontal slice and vertical cross-section of # almost uniformly polarized structure and striped domain structure. #== Proposal of an experiment of BaTiO3 thin-film capacitor Fig. 1 figures/proposal2.gif (1) On SrRuO3 or Pt metallic substrate, (2) grow BaTiO3 thin and thick thin films, (3) cover them with the metal, (4) if you use x-ray or whatever, (5) you may observe deference of fineness in their domain structures. /Fig. == Applied external electric field Fig. 1 figures/step.jpg Schematic illustrations of triangle-wave electric field used to measure ferroelectric hysteresis loops experimentally (inset) and in the present simulations. /Fig. Fig. 1 figures/external-E-field-hysteresis-loop.jpg Latest version of feram can apply smooth triangle-wave electric field. /Fig. == Temperature dependence of hysteresis loops for bulk BaTiO3 Fig. 1 figures/bulk-two-box.jpg Temperature dependence of hysteresis loops for bulk BaTiO3. The phase transition from paraelectric phase to ferroelectric phase is clearly seen. /Fig. == Hysteresis loops for epitaxially constrained and "free" BaTiO3 film capacitors Fig. 1 figures/epit-vs-free-hysteresis-box.jpg Calculated hysteresis loops for capacitors with (a) epitaxially constrained films, and (b) ``free'' films of various thickness $l$ and with dead layer $d$ at 100 K. /Fig. # Supercell sizes were $16\times 16\times 2(l+d)$. # $\Delta t=2$~fs and $n_{\rm steps}=$~50,000. # For epitaxially constrained films, # $\mathcal{E}_0 =$~4,000~kV/cm and $\Delta \mathcal{E} =$~100~kV/cm are employed. # For ``free'' films, # $\mathcal{E}_0 =$~400~kV/cm and $\Delta \mathcal{E} =$~10~kV/cm are employed. == Out-of-plane polarization is no longer the ground state Fig. 1 figures/damaged.jpg Vertical cross section of a simulated ferroelectric "free" film capacitor with a single dead layer; 16 x 16 x (l=63, d=1). The snapshot was taken at the point marked "x". The projection of the dipole moments onto the xz-plane are indicated with arrows. /Fig. == Epitaxially constrained film vs. "free" film (again) Fig. 1 figures/epit-vs-free.jpg Schematic comparison between the epitaxially constrained film and the ``free'' film. In the epitaxially constrained film, switching may have to climb over a potential barrier, but, in the ``free'' film, dipoles can be easily rotated and switching can go around a valley of the potential. /Fig. == Future plan with the feram code: Molecular dynamics simulations of relaxors Relaxors Pb$B'_x$$B''_{1-x}$O$_3$: Pb$^{2+}$ * Pb(Sc$_{1/2}$Nb$_{1/2}$)O$_3$ (PSN): Sc$^{3+}$, Nb$^{5+}$ * Pb(Sc$_{1/2}$Ta$_{1/2}$)O$_3$ (PST): Sc$^{3+}$, Ta$^{5+}$ * Pb(Mg$_{1/3}$Nb$_{2/3}$)O$_3$ (PMN:) Mg$^{2+}$, Nb$^{5+}$ * $(1-x)$Pb(Mg$_{1/3}$Nb$_{2/3}$)O$_3$-xPbTiO$_3$ (PMN-$x$PT): Mg$^{2+}$, Nb$^{5+}$, Ti$^{4+}$ * Pb(Mg$_{1/3}$Ta$_{2/3}$)O$_3$ (PMT): Mg$^{2+}$, Ta$^{5+}$ * Pb(Zn$_{1/3}$Nb$_{2/3}$)O$_3$ (PZN): Zn$^{2+}$, Nb$^{5+}$ The averaged valence number of B-site ions $B'_x$ and $B''_{1-x}$ is $4+$. # There are some Pb-based relaxors PbB'B''O$_3$. == Frequency dependence of dielectric constant of relaxors Fig. 1 figures/TsurumiFig3.jpg After [Takaaki Tsurumi, Kouji Soejima, Toshio Kamiya and Masaki Daimon: Jpn. J. Appl. Phys. 33, 1959-1964 (1994)] /Fig. # Most characteristic and still mysterious property of relaxors is #Eq.1 #\tan\delta=\tan() #/Eq. == Crystal structures of relaxors Fig. 1 figures/RandallBhallaFig3.jpg (a) Ordered B-site structure of $AB'_xB''_{(1-x)}$O$_3$. In the case of Pb(Sc$_{1/2}$Nb$_{1/2}$)O$_3$, Sc and Nb has the NaCl structure. In the case of Pb(Mg$_{1/3}$Nb$_{2/3}$)O$_3$, the Na-site is randomly occupied by Mg$_{1/3}$Nb$_{1/6}$ and the Cl-site is occupied by Nb$_{1/2}$. (b) Disordered structure. B-sites is occupied by $B'$ and $B''_{(1-x)}$ randomly. After [C. A. Randall and A. S. Bhalla: Jpn. J. Appl. Phys. 29, 327-333 (1990)]. /Fig. == Local field on Pb-site Displacement of Pb is the main source of dipole moment. Fig. 1 figures/TwoDoubleWells.jpg Symmetric potential around Pb at the chemically ordered region (COR). Asymmetric one at the chemically disordered region (CDR). /Fig. == Chemically ordered and disordered regions of relaxors Fig. 1 figures/RandallBhallaFig8.jpg The quenched sample has chemically ordered regions (∼100Å) in a chemically disordered matrix. After [C. A. Randall and A. S. Bhalla: Jpn. J. Appl. Phys. 29, 327-333 (1990)]. /Fig. == Modeling of relaxors Fig. 1 figures/DoubleWells.jpg Schematic illustration of modeling of relaxor. Bubble-like Chemically ordered regions (COR, red symmetric double wells) is created in a chemically disordered matrix (blue asymmetric double wells) /Fig. == Snapshot of a simulation of PSN relaxor Fig. 1 figures/SilviaFig1.jpg A (110) plane through the PSN simulation box representing the projected local field (arbitrary units) at each Pb site in the plane. Chemically ordered regions (approximately circular) have small approximately homogeneous fields, and chemically disordered regions have larger more varied and disordered local fields. After [Tinte, Burton, Cockayne and Waghmare: Phys. Rev. Lett. 97, 137601 (2006)]. /Fig. == Frequency dependence of dielectric constant Eq. (Fluctuation-dissipation theorem) \epsilon(\omega) - \epsilon_\infty = \frac{1}{3\epsilon_0 V k_B T}\left[ \left\langle \bm{p}^2 \right\rangle + i\omega \int_{0}^{\infty} dt e^{i\omega t} \left\langle \bm{p}(t) \cdot \bm{p}(0) \right\rangle \right] /Eq. where $\bm{p}(t)$ is the total electric dipole moment in the supercell at time $t$, Eq. (Total electric dipole moment in the supercell) \bm{p}(t)=Z^*\sum_{\bm{R}} \bm{u}(\bm{R};t) . /Eq. Estimation of computational time 32x32x32 unit cells, $\Delta t = 2$ fs, [AMD64 1.8GHz dual core] or [SR11000 1 node = 16 cores] 1THz, $T=1\,$ps, 500 steps, 54 sec or 8.4 sec. 1GHz, $T=1\,$ns, 500,000 steps, 900 min. or 140 min. 1MHz, $T=1\mu$s, 500,000,000 steps, 620 days or 97 days So, I am planning to calculate and compare $\epsilon(\omega)$ for $\omega/(2\pi)$ = 10M ... 10GHz. # If you apply external electric field at site R' at time t', dipole at R,t responses # through this equation with real-space-and-time susceptibility tensor chi(R-R',t-t'). # You can simplify it by Fourier transformation and get q-omega expression. # And then, if you apply the fluctuation-dissipation theorem to it, you may get # the imaginary part of chi. Then with Kramers-Kronig relation, you can get the real part. # # G. L. Squires: Introduction to the Theory of Thermal Neutron Scattering (Dover Publications, 1996). == Summary * We developed "feram", a fast simulator for perovskite-type ferroelectric bulks and thin films * Molecular dynamics (MD) simulation with first-principles-based effective Hamiltonian * Phase transition of bulk ferroelectrics * Thin-film capacitor: perfect and imperfect short-circuited electrodes * Striped domain structure is predicted * Hysteresis loops for epitaxially constrained and "free" BaTiO3 film capacitors * Investigations of relaxors with MD are proposed. Fig. 1 ../logo.jpg Thank you for your attention. /Fig. # * An experiment of BaTiO3 thin-film capacitor is proposed =end # # for ulmul.rb # require 'rubygems' require 'ulmul' u=Ulmul.new(1..0,'ulmul2xhtml') u.subs_rules << [ /(heating-up)\s*(\S*\/?\.gif)(\s|$)/,'\1\3'] u.subs_rules << [/(cooling-down)\s*(\S*\/?\.gif)(\s|$)/,'\1\3'] u.parse(ARGF) puts u.html(["feram-presentation.css"],["ulmul-slidy.js"],"Takeshi Nishimatsu","en") # Local variables: # mode: RD # compile-command: "ruby film.txt film.txt > film.xhtml" # End: