Oxide Dissolution Mediated Formation of Conducting Filament in ReRAM Devices: A Phase Field Study

Abstract

Formation of oxygen vacancy-rich conducting filament in the metal-oxide insulating layer due to applied electric field is studied using phase field modelling. We consider that the formation of the Frenkel defects plays a dominant role in the process of filament formation. The choice of any other type of defects have been neglected due to the consideration of charge neutrality. We find that during the initial stages of filament formation, intrinsic defect concentration plays a crucial role. However, during the late stages, when the filament formation approaches completion, because of the increased electric field, the generation of new Frenkel defects dominates the growth mechanism. Our observation confirms that the intrinsic defects are not sufficient for the completion of conducting filament in the insulating memresistive layer. Our numerical analysis helps us to better understand the formation mechanism of conducting filament and to determine the key material parameters that influence the operation of oxide based non-volatile random access memory (RAM) i.e. resistive RAM (ReRAM).

1. Introduction

In recent years, the need for a new generation of non-volatile random access memory (RAM) increases the academic and the industrial attentions for the resistive switching based RAM (ReRAM) devices^1–9) substantially. The key mechanism that drives the electric responses in the ReRAM devices is the formation of the conducting filament (CF). Also, the way CF forms during the operation conditions, vary with the choice of memresistive systems. For example, migration of positively charged cation could result in the formation of metallic CF¹⁰⁾ or negatively charged anion migration could result in the formation of metal-oxide CF.¹¹⁾

Moreover in the case of a few metal-oxide systems, such as ZnO, the mechanism of the formation of CF is different. In such a case, according to the previous experimental work of Chen et al.,¹⁰⁾ oxygen vacancy plays an important role to form the Zn rich CF. Intrinsic defect concentration, namely oxygen vacancy, is responsible for the formation mechanism.

A number of continuum models have been employed to predict the morphological changes that occur in the electrochemical system. These computational tools are based on the phase field modelling and the drift diffusion modelling. Using the phase field models developed by Guyer et al.^12,13) and Shibuta et al.^14,15) to study electrochemical processes, we could study the switching process in detail. We have implemented numerical techniques to understand the related electrical response in model ReRAM devices.^16–18) Formation of CF consisting of positively charged species in a Cu-based switching layer has been developed by Shibuta et al.¹⁵⁾ and Okajima et al.¹⁹⁾ Typical I-V behaviour due to voltage sweeping has been extracted successfully from the model. It is important to mention that these phase field models are also well suited to study the electro deposition process. For this reason, the incoming boundary flux of charged species is an essential part of these models. Phase field model to study the morphology of deposited Li in Li-based batteries have been proposed by Liang et al.^20,21) Here, the deposited Li also develops the filamentary shape during the underlying electrochemical processes. Drift-diffusion models are also employed to study the formation of CF in oxide based memresistive devices.²²⁾ These models are particularly useful in predicting the I-V behaviour during various operation conditions.

We are interested in investigating the role of electric field driven phase transformation in the formation mechanism. In metal-oxide insulating layer, such as CuO²³⁾ and NiO,²⁴⁾ the cationic metal ion or vacancy contributes in forming the CF. However, oxygen vacancy also plays a dominant role in developing the CF under applied electrical bias in various metal-oxide systems. The memresistive properties of Cr doped SrTiO3 are observed due to the creation of Cr-rich metallic filament formed due to the electric field induced migration of oxygen vacancies.²⁵⁾ Oxygen vacancies are also crucial in obtaining switching behaviour in TiO₂ based ReRAM devices.^26,27) On the other hand in Ta₂O_5−x/TaO_2−x, oxygen vacancies are typically supplied from the TaO_2−x layer and produce switching behavior by forming a conducting path.²⁸⁾ Using transmission electron microscope (TEM), Chen et al.¹⁰⁾ have observed the formation of Zn-rich CF in the ZnO system during the voltage sweeping process. They suggest that due to the electric field driven phase transformation, dissolution of ZnO occurs, which produces positively charged oxygen vacancies. Due to repulsion of positively charged defects towards the low electric field, they accumulate at the top electrode and produce Zn-rich CF. Using the phase field model we can consider such a system where CF is composed of a metal-rich phase and a metal-oxide matrix consisting of intrinsic oxygen vacancies.

However, as we have discussed in this paper, oxygen vacancies may not be the only key player. We consider that the other charged species, which are obvious for this type of systems, such as Frenkel pairs of negatively charged oxygen interstitials and already considered oxygen vacancies, are coexisting in the system. Also, due to the applied electric field new vacancies and interstitials can be generated using the bond-breaking model.^29,30) The rate of creation of new defects can be described in terms of the activation energy and the defect formation energy. The role of such materials parameters prove to be crucial in fine tuning the device applicability of memristors. In this study, we explore the role of such materials parameters in detail.

The contents of this paper are as follows: we start with the description of the phase field model, used in exploring the formation of CF. Detailed description of the phase field, concentration of positively charged vacancy and negatively charged interstitial have been described. The rate of formation and annihilation of Frenkel defects have been incorporated in the governing equations for defect concentrations. In the result section, we explain the electric field driven growth or dissolution of the metal-rich CF in the metal-oxide insulating layer. The role of intrinsic defect concentration and its effect on the growth of CF are studied. From our numerical results, we could justify the need to explore the role of Frenkel defect formation due to the applied electric field. After incorporating the growth morphology supported by the Frenkel defect and applied electric field, we investigate the role of material parameters on typical I-V behaviour of memristor devices. We also investigate the role of operation condition, such as the effect of voltage sweeping rate on the I-V curve. It is important to mention that the proposed electric field driven growth mechanism is general for such memristor systems where the dissolution of switching layer and the creation of defects play a dominant role. Hence, we can employ this model for a wide variety of ReRAM systems and explore the effect of different material and operation parameters on the device performance.

2. Phase Field Modelling

We consider that the resistive switching layer is composed of a metal-oxide matrix; consisting of positively charged intrinsic oxygen vacancy (V⁺). We denote this phase as the α matrix and the concentration of oxygen vacancy is described by $C_{V^{ + }}^{\alpha }$. Also, we consider that the matrix phase is supersaturated with V⁺. The system with such a supersaturated matrix increases the driving force for nucleation and also helps in nucleating a vacancy-rich β-phase precipitate in the metal-oxide matrix. We denote the concentration of the charge neutral vacancy-rich phase as $C_{M}^{\beta }$. Using ideal solution model, we write the bulk free energies of α-phase (f^α) and β-phase (f^β) as,   

\begin{align} f^{\alpha} & = \left(\text{G}_{B}^{\alpha} + \frac{RT}{V_{m}}\ln\frac{1 - {}^{e}C_{V^{+}}^{\alpha}}{{}^{e}C_{V^{+}}^{\alpha}} \right) C_{V^{+}}^{\alpha} + \text{G}_{B}^{\alpha} (1 - C_{V^{+}}^{\alpha})\\ & \quad + \frac{RT}{V_{m}}(C_{V^{+}}^{\alpha} \ln C_{V^{+}}^{\alpha} + (1 - C_{V^{+}}^{\alpha}) \ln (1 - C_{V^{+}}^{\alpha})) \end{align} (1)

  

\begin{align} f^{\beta} & = \left(\text{G}_{B}^{\beta} + \frac{RT}{V_{m}}\ln\frac{1- {}^{e}C_{M}^{\beta}}{{}^{e}C_{M}^{\beta}} \right) C_{M}^{\beta} + \text{G}_{B}^{\beta} (1 - C_{M}^{\beta})\\ &\quad + \frac{RT}{V_{m}} (C_{M}^{\beta} \ln C_{M}^{\beta} + (1 - C_{M}^{\beta}) \ln (1 - C_{M}^{\beta})) \end{align} (2)

Here, ${}^{e}C_{V^{ + }}^{\alpha }$ (0.5) and ${}^{e}C_{M}^{\beta }$ (0.98) represent the equilibrium concentration of α and β phases, respectively. R is the gas constant. T = 300 K and V_m = 1.4517 × 10⁻⁵ m³ mol⁻¹ are the temperature and the molar volume of the system. $G_{B}^{\alpha } = 5.199$ RT/V_m and $G_{B}^{\beta } = 1.98$ RT/V_m describe the thermodynamic parameters for α and β phases, respectively. For a particular supersaturation of $C_{V^{ + }}^{\infty }$, we can calculate the driving force for growth. With this assumption – unless the additional defects are generated in the system – a β vacancy-rich precipitate of ≈2.0 nm of (critical) radius grows as long as the $C_{V^{ + }}^{\infty }$ (if chosen to be 0.51) attains the equilibrium matrix concentration of ${}^{e}C_{V^{ + }}^{\alpha }$, see our previous publication¹⁷⁾ for an elaborate discussion. Moreover, the applied electric potential also modulates the dynamics of the positively charged V⁺ in the system. Due to such an application of electric potential, positively charged V⁺ migrates toward the low electric field. Hence it is easy to comprehend the fact that the supersaturation of the matrix phase (which is intrinsically supersaturated with the V⁺) combined with the applied electric potential govern the growth of the filament/matrix interface. The total free energy density of such a system can be expressed as,   

\begin{equation} f = \int_{V} \left( f_{\textit{bulk}} + \frac{\kappa_{\phi}{}^{2}}{2} |\nabla \phi|^{2} + \frac{1}{2}\rho \Phi \right)dV \end{equation} (3)

Here ϕ describes the phase field – to distinguish the α and the β phases – of the system. Specifically, ϕ attains a value of 0 in the α phase and a value of 1 in the β phase. The gradient free energy coefficient is denoted by κ_ϕ ($ = \sqrt{\frac{6\lambda \sigma }{\alpha_{\lambda }}} $). σ = 0.025 J m⁻² and 2λ = 1.0 × 10⁻⁹ nm are the model parameters representing the interfacial energy and the interfacial width of the system. α_λ = 2.94 is a parameter which depends on the choice of range in ϕ to distinguish the interfacial region and is typical in the Kim-Kim-Suzuki (KKS) phase field model.^31,32) The bulk free energy and the electrostatic energy of the system are expressed by f_bulk and $\frac{1}{2}\rho \varPhi $ respectively. Here, electric charge of the species is denoted by $\rho = \frac{Fz_{V^{ + }}}{V_{m}}C_{V^{ + }}^{\alpha }$ (F is the Faraday constant and $z_{V^{ + }} = + 2$ is the valence of the V⁺) and the applied electric potential of the system is denoted by Φ. By following the KKS phase field model^31,32) for binary alloy and using the bulk free energy of α (eq. (1)) and β (eq. (2)) phases, f_bulk can be written as,   

\begin{equation} f_{\textit{bulk}} = (1 - p (\phi))f^{\alpha} (C_{V^{+}}^{\alpha}) + p (\phi) f^{\beta} (C_{M}^{\beta}) + Wg(\phi) \end{equation} (4)

Here p(ϕ) = ϕ³(10 − 15ϕ + 6ϕ²) is known as the Wang interpolation function,³³⁾ which attains a value of 0 in the α-matrix (ϕ = 0) and 1 in the β-precipitate (ϕ = 1) with a continuous variation in the 0 ≤ ϕ ≤ 1 range. The energy barrier between matrix and precipitate phases are set by the constant $W = \frac{3\sigma \alpha_{\lambda }}{\lambda }$ J m⁻³ times the double well potential, g(ϕ) = ϕ² (1 − ϕ)². We provide the detailed description of the terms used so far in our previous publication.¹⁷⁾ In Refs. 17), 18), we consider the filament formation in the presence of supplied charged species from the active electrode as has been reported in case of electrochemical metallization (ECM) memory cells. We also assume that the growth of β precipitate is supported by the formation and annihilation of vacancy and interstitial Frenkel defect pair in the system. In this case, the formation of defects is assisted by the applied electric field. The rate can be described by,²⁹⁾   

\begin{equation} k_{f} = A_{f}\exp \left(-\frac{E_{a} + E_{f} - bE_{\Phi}}{RT}\right) \end{equation} (5)

A_f is a constant (set to unity in numerical simulations) and it describes the rate factor. The activation energy and the defect formation energy are presented in terms of E_a and E_f respectively. Due to the applied electric field E_Φ, the energy required to break the di-electric bond of metal-oxide is bE_Φ; b is known as the bond-polarization factor.^29,30) Similarly, the rate of annihilation of Frenkel defects can be expressed using,   

\begin{equation} k_{a} = A_{a}\exp \left(-\frac{E_{a}}{RT} \right) \end{equation} (6)

Here, A_a is the annihilation rate factor and is set to unity in numerical simulations similar to A_f. We schematically describe the mechanism of initial filament formation and subsequent creation of Frenkel defects due to applied electric potential in the system to support the growth of CF in Fig. 1. The temporal evolution of phase field can be obtained from the chemical potential; deduced by calculating the variation of total free energy density f (eq. (3)) with respect to ϕ. The governing equation for ϕ can be written as,   

\begin{align} \frac{\partial \phi}{\partial t} &= M_{\phi}\bigg[\kappa_{\phi}{}^{2}\nabla^{2}\phi \\ &\quad+ \frac{dp(\phi)}{d\phi}\frac{RT}{V_{m}}\left(f_{C_{V^{+}}{}^{\alpha}}^{\alpha} - f_{C_{M}{}^{\beta}}^{\beta} - \frac{\partial f_{C_{V^{+}}}^{\alpha}}{\partial C_{V^{+}}{}^{\alpha}} (C_{V^{+}}{}^{\alpha} - C_{M}{}^{\beta}) \right) \\ &\quad- W\frac{dg(\phi)}{d\phi} \bigg] \end{align} (7)

Fig. 1

Variation in Gibbs free energy as a function of vacancy composition is shown in (a). The initial oxygen vacancy composition, intrinsic in nature, provides the necessary driving force for vacancy-rich conducting filament (CF) phase to nucleate at the bottom electrode and is shown in (d). The available supersaturation or the intrinsic vacancy composition drive the filament growth under the influence of positive bias applied in the top electrode. As the filament grows, the separation distance between tip of the CF and the top electrode reduces which in turn produce the intensification of electric field, as shown schematically in (b). Such intensification in the electric field promotes the dielectric breakdown of the oxide phase and generates the Frenkel pair of oxygen vacancies and interstitials. The energetics associated with overcoming the activation energy (E_a) and modulating the defect formation energy (E_f) through the bond-polarization (b) times the electric field (E_Φ) is presented schematically in (c).

Considering the formation and annihilation of vacancy and interstitial, the temporal evolution can be equated with the diffusional flux as follows:   

\begin{align} \frac{\partial C_{V^{+}}}{\partial t} & = \nabla \left[\frac{D (\phi)}{f_{cc}}\nabla \left(f_{C_{V^{+}}}^{\alpha} + \frac{z_{V^{+}}F\Phi}{2} \right) \right]\\ &\quad + k_{f} (1 - C_{V^{+}}) (1 - C_{I^{-}}) - k_{a} C_{V^{+}}C_{I^{-}} \end{align} (8)

  

\begin{align} \frac{\partial C_{I^{-}}}{\partial t} & = \nabla \left[\frac{D_{I^{-}}}{f_{cc}{}^{I^{-}}}\nabla \left(f_{C_{I^{-}}} + \frac{z_{I^{-}}F\Phi}{2} \right) \right]\\ &\quad + k_{f}(1 - C_{V^{+}}) (1 - C_{I^{-}}) - k_{a} C_{V^{+}} C_{I^{-}} \end{align} (9)

The phase field mobility and diffusivity of the oxygen vacancy are termed as M_ϕ and D(ϕ) (= D_βp(ϕ) + D_α(1 − p(ϕ))), respectively. The value for phase field mobility is set to 7.7424 × 10⁻⁸ mol J⁻¹s⁻¹. We choose the diffusivity of oxygen vacancy in α and β phases to be D_α = 3.0 × 10⁻¹⁸ m²s⁻¹ and D_β = 3.0 × 10⁻²⁰ m²s⁻¹ respectively. We assume the diffusivity of the interstitial to be $D_{I^{ - }} = 3.0 \times 10^{ - 20}$ m²s⁻¹ in the matrix and in the conducting phase, we assume the diffusivity to be zero. The electric potential distribution in the system can be deduced by solving the following Poisson’s equation of charge conservation,   

\begin{equation} \nabla \cdot [\varepsilon (\phi) \nabla \Phi] = 0 \end{equation} (10)

Here, ε(ϕ) (= ε_βp(ϕ) + ε_α(1 − p(ϕ)) is the conductivity of the system. ε_α = 2.0 S m⁻¹ and ε_β = 2.0 × 10⁴ S m⁻¹ are the conductivities of the α and the β phases respectively.

3. Numerical Implementation

Second order accurate finite difference scheme has been implemented to solve the governing equations for phase field, vacancy and interstitial concentrations. In order to obtain the electric potential distribution using eq. (10), we employ the Gauss-Seidel method. To obtain the value of diffusivity in case of eq. (8) and eq. (9), and conductivity in case of eq. (10) at the interface (0 ≤ ϕ ≤ 1), we employ the sixth order Lagrange interpolation formula. We choose the simulation domain in such a way that the left and the right boundary represent the counter electrode and active electrode, respectively. For this reason, we adopt the adiabatic boundary condition along these two boundaries. On the other hand, along the top and the bottom boundary periodic boundary conditions are applied. Hence, the simulation domain can be considered to be a long channel of resistive layer sandwiched between two electrodes. The numerical code is implemented using CUDA-C for GPU parallelization on Nvidia Titan Xp graphics.

4. Results and Discussions

In order to better understand the dynamics of the charged particle with the application of electric potential, we consider a planar interface located at the middle of the simulation domain. The conducting phase and the matrix phase consisting of oxygen vacancies are located at the left and the right of the planer interface, respectively. The composition of the conducting phase is set to 0.98 and the matrix phase consisting of oxygen vacancy is set to 0.5, which by the choice of Gibbs free energy represent the equilibrium composition for the corresponding phases. The left boundary of the simulation domain represents the bottom/counter electrode and the right boundary of the simulation domain represents the top/active electrode. Here, we neglect the formation and annihilation terms which control the generation of the Frenkel defects in the system. With the application of electric potential of positive polarity on the right boundary (or top electrode), positively charged vacancies move towards the region of low electric potential (the bottom electrode is kept at 0 V) and the interface of the conducting phase advances towards the right boundary – see Fig. 2(a) for cases with +ve applied voltage. On the other hand, with the application of negative electric potential on the right boundary, dissolution of the metal-rich conducting phase is observed and the interface shrinks away from the right boundary – see Fig. 2(a) for cases with −ve applied voltage. We can see from Fig. 2(a), that in case of Φ = 0 V, the concentration at the right boundary remains at its equilibrium value (which is 0.5). It happened, as we have set the far field concentration of the system to its equilibrium value while examining the effect of applied voltages. For this reason, because of the absence of driving force (thermodynamics and electric field driven), interface migration can not be observed. However, with varying polarity and strength of the electric field, it attains a specific value above (due to negative electric potential) or below (due to positive electric potential) the equilibrium compositions at the active electrode or the simulation boundary at the right. From Fig. 2(b), we understand that the amount of dissolution or accumulation of vacancy concentration at the right boundary depends on the polarity and strength of the applied electric potential. Nevertheless, the migration of metal rich precipitate phase stops as soon as the vacancy concentration on the right boundary reaches a specific value, which is determined by the polarity and magnitude of the applied electric voltage – see Fig. 3(a) and (b) for the composition and electric potential profile across the interfaces due to −ve and +ve applied voltages, respectively. At this point, we want to note that the obtained composition distributions for the planar interface profile are consistent with the numerical results described by Shibuta et al.^14,15)

Fig. 2

(a) Variation of oxygen vacancy composition across a planer interface obtained for different choices of applied electric voltage. Attained composition at the right boundary (top electrode) due to the applied electric voltages are shown in (b). The figures are produced by neglecting the additional Frenkel defect formation and annihilation terms.

Fig. 3

Oxygen vacancy composition distribution across a planer interface obtained by neglecting the additional Frenkel defect formation and annihilation terms with application of negative voltage is shown in (a) and in case of applied positive voltage is shown in (b).

Hence, it is easy to comprehend the fact that the system with only intrinsic vacancy concentration ceases to grow after it depletes in the system. With this observation, we came to the conclusion that intrinsic vacancy concentration is not sufficient to complete the formation of the filament. Due to such reasons, we choose the growth model with the formation and annihilation of the Frenkel defect due to the applied electric field. However, we can see from eq. (5) that due to application of a very large electric field, the exponential term may attain undefined numerical values. Hence, we stop our simulation at a moment when the electric field at the α − β boundary exceeds 1.5 × 10⁸ V/m. We present the formation of the conducting filament in Fig. 4.

Fig. 4

Electric potential induced morphological changes of the conducting filament and the corresponding I-V curve are shown in the figure. The distribution of Frenkel defect pairs of oxygen vacancy and interstitials along with the distribution of electric potential are presented. Intrinsic oxygen vacancy assists the growth of CF till point (2). For further growth of the filament, additional Frenkel defect pairs are generated with the increasing magnitude of electric field in the system and assist the CF to grow further. At point (4), due to the generation of oxygen interstitials, the system depletes in oxygen vacancy and limits the further growth of the CF.

We study the morphological evolution of CF to understand the parameters that influence the switching phenomena in oxide based resistive RAM. The corresponding I-V responses – due to the application of different rates of voltage sweeping – are recorded during the growth of the filament. The voltage at which the current value increases abruptly is known as the forming voltage. In this paper we explore the effect of material parameters as well as operation conditions on the device performance in terms of such recorded forming voltages. A complete switching process in ReRAM devices, consisting of forming, SET, and RESET stages will be discussed in studies which will be published at some later time. We start our discussion with the effect of the bond-polarization factor on the forming voltage.

In our model, we have assumed that the electric field generates new Frenkel defects during the operation of the oxide based ReRAM devices. As a result, the generated oxygen vacancies contribute to the formation of CF and promote the electrical responses. To better understand the effect of variation in such parameters which controls the defect concentration, we calculate the I-V curve under various operation conditions and have studied the changes in the forming voltages. We choose three different numerical values for b; 30.1 eA, 40.4 eA and 50.8 eA. I-V curves for these three different choices of b are shown in Fig. 5(a). The estimated current values do not deviate for these three cases till ≈0.32 V. This also indicates that the growth rate is the same at the initial stages for three different choices of b. However, for the voltages >0.32 V, the measured current value deviates considerably. In particular, it is observed that the forming voltage increases with the decreasing value of b, see Fig. 5(b). Such variation in the forming voltage manifests from the fact that the reduced magnitude of b limits the generation of new Frenkel defects (prerequisite for the growth of CF) considerably, which in turn increases the forming voltage. At low operating voltages, due to the difficulties in generating the new defects, the growth of filaments are mainly controlled by the intrinsic defects (which are oxygen vacancy in this case) of the system. Because of such processes, at low operating voltages, the growth of CF is similar for the different choices of b. We term this region as “intrinsic defect-assisted”. On the other hand, at relatively high operating voltages, the growth of CF is dominated by the newly generated Frenkel defects, which is termed as “Frenkel defect-assisted”.

Fig. 5

Effect of bond-polarization on I-V curve is shown in (a). The variation of forming voltage with the choice of bond-polarization is shown in (b).

Formation of new charged species – crucial for the growth of CF – depends on the choice of the memresistive system, in which case the choice of activation energy (E_a) and defect formation energy (E_f) may vary significantly. Using our phase field model, we can indeed study how filament forming voltage varies with such choices, which control the rate of defect formation under applied electric stimulus. We show the I-V curves for different choices of E_f in Fig. 6(a) and in Fig. 6(b), we show how the forming voltage changes with the choice of E_f. Here, we note that the E_a is set to 0.106 eV and we record the I-V curve for different choices of E_f, as shown in Fig. 6. The energy barrier – which is required to overcome in order to form sufficient defects – is reduced in cases of low magnitude of E_f. Due to such favourable conditions of defect formation, newly formed defects support the filament growth and reduce the forming voltage when E_f is relatively low.

Fig. 6

Effect of defect formation energy on I-V behaviour is presented in (a). In (b), we show the variation of forming voltage with changing formation energy.

We record the I-V curve for different choices of activation energy (E_a) and we set the formation energy (E_f) to 0.149 eV in order to explore the effect of E_a on the forming voltage. We find from Fig. 7(a) that the forming voltage increases substantially with the increase of E_a till a specific value (as set by the choice of b) after which the increase in forming voltage is reduced considerably. This could be attributed to the fact that the choice of E_a sets the limit on the filament growth kinetics. An increased value of E_a, increases the activation barrier of defect formation and thereby reduces the filament growth. In Fig. 7(b), we show the defect formation rate (k_f) as has been described in eq. (5) at the filament tip when the formation of CF completes for a varying choice of E_a. k_f significantly reduces with the increasing magnitude of E_a. We can see from Fig. 7 that k_f becomes negligible as the magnitude of E_a increases and thereby limiting the supply of defects crucial for the growth of CF. For this reason, we understand that a metal-oxide memresistive system with relatively high magnitude of activation energy for defect formation is not desirable from the perspective of the device operation. Such a system may impose limitation in the growth of the vacancy-rich conducting filament and thereby limits the desired switching characteristics of the devices.

Fig. 7

In (a), we show how the choice of activation energy modulates the forming voltage in systems with three different values of bond-polarization. The defect formation rate at the tip of the filament for different choices of activation energy is presented in (b).

We have also explored the effect of operation conditions on the forming voltage. We choose a varying rate of voltage sweeping as the operation parameter. We consider three different rates of voltage sweeping; 0.01 V/s, 0.05 V/s and 0.1 V/s. The resulting I-V curve is shown in Fig. 8. We find that the forming voltage increases with the increasing rate of voltage sweep. Also, the ratio of forming voltage and sweeping rate describes the forming time of the CF under applied electrical stimulus. In this case, we can understand that the increase in sweeping rate actually indicates the decrease in filament forming time. It is important to note that the filament growth is not only controlled by the availability of the defects in the system – as have been discussed in the previous sections – but also by the migration of it to the tip of the growing filament. Due to such a kinetic process, with the increasing magnitude of applied electrical bias (caused by the increasing voltage sweeping rate), kinetics of charge migration dominate and forming time decreases with the increasing rate of sweeping voltage.

Fig. 8

Effect of voltage sweep rate on I-V curve.

5. Conclusion

We implement the mechanism of CF growth as has been described by Chen et al.¹⁰⁾ using our phase field model. From the morphologies presented in Fig. 4 and the observation presented in Fig. 5, we understand that only intrinsic defect, which is oxygen vacancies in this case, may not be sufficient to complete the formation of CF and for this reason we have incorporated a mechanism to generate additional Frenkel defects to support the filament growth. However, in such cases where defects are generated in pairs, it becomes challenging to overcome the separation distance between the top electrode and the growing CF at the late stages, due to the enrichment of the oxygen interstitials at high applied electric field. Such an enriched region of oxygen interstitial depletes the oxygen vacancies and hinders the growth of filament. We believe that such challenges in describing the completions of filament formation opens up the possibility of new postulation in conduction mechanism when the separation distance is reduced considerably to ≈2.0 nm. In such a situation, the system cannot provide sufficient oxygen vacancies to complete the formation of CF. Moreover, even with the challenges we faced in obtaining the completely formed CF which essentially short-circuits the top and bottom electrodes, we explore the effect of b, E_a, and E_f on the electric responses in the model system. We connect the filament formation mechanism in terms of morphological evolution of CF with the energetics (described using E_a and E_f of the system) and electric field induced modulation (bE_Φ) of defect formation.

Acknowledgments

We thank the funding agencies for providing essential support to conduct this research. We acknowledge Creative Materials Discovery Program through the National Research Foundation of Korea (NRF) funded by the Ministry of Science, ICT and Future Planning (2016M3D1A1027666), and National Research Foundation of Korea (NRF) grant funded by the Korea government (MSIP) (NRF - 2019R1A2C1089593 and NRF - 2019M3A7B9072144) for funding support.

REFERENCES

• 1)   R.  Waser,  R.  Dittmann,  C.  Staikov and  K.  Szot: Adv. Mater. 21 (2009) 2632–2663. doi:10.1002/adma.200900375
• 2)   L.  Zhu,  J.  Zhou,  Z.  Guo and  Z.  Sun: J. Materiomics 1 (2015) 285–295. doi:10.1016/j.jmat.2015.07.009
• 3)   D.  Ielmini: Semicond. Sci. Technol. 31 (2016) 063002. doi:10.1088/0268-1242/31/6/063002
• 4)   M.A.  Zidan,  A.  Chen,  G.  Indiveri and  W.D.  Lu: J. Electroceram. 39 (2017) 4–20. doi:10.1007/s10832-017-0103-0
• 5)   Z.  Zhang,  Z.  Wang,  T.  Shi,  C.  Bi,  F.  Rao,  Y.  Cai,  Q.  Liu,  H.  Wu and  P.  Zhou: InfoMat 2 (2020) 261–290. doi:10.1002/inf2.12077
• 6)   Y.  Chen: IEEE Trans. Electron Devices 67 (2020) 1420–1433. doi:10.1109/TED.2019.2961505
• 7)   T.  Li,  H.  Yu,  S.H.Y.  Chen,  Y.  Zhou and  S.T.  Han: J. Mater. Chem. C 8 (2020) 16295–16317. doi:10.1039/D0TC03639K
• 8)   F.  Zahoor,  T.Z.  Azni Zulkifli and  F.A.  Khanday: Nanoscale Res. Lett. 15 (2020) 90. doi:10.1186/s11671-020-03299-9
• 9)   E.  Carlos,  R.  Branquinho,  R.  Martins,  A.  Kiazadeh and  E.  Fortunato: Adv. Mater. 33 (2021) 2004328. doi:10.1002/adma.202004328
• 10)   J.Y.  Chen,  C.L.  Hsin,  C.W.  Huang,  C.H.  Chiu,  Y.T.  Huang,  S.J.  Lin,  W.W.  Wu and  L.J.  Chen: Nano Lett. 13 (2013) 3671–3677. doi:10.1021/nl4015638
• 11)   J.  Shang,  G.  Liu,  H.  Yang,  X.  Zhu,  X.  Chen,  H.  Tan,  B.  Hu,  L.  Pan,  W.  Xue and  R.W.  Li: Adv. Funct. Mater. 24 (2014) 2171–2179. doi:10.1002/adfm.201303274
• 12)   J.E.  Guyer,  W.J.  Boettinger,  J.A.  Warren and  G.B.  McFadden: Phys. Rev. E 69 (2004) 021603. doi:10.1103/PhysRevE.69.021603
• 13)   J.E.  Guyer,  W.J.  Boettinger,  J.A.  Warren and  G.B.  McFadden: Phys. Rev. E 69 (2004) 021604. doi:10.1103/PhysRevE.69.021604
• 14)   Y.  Shibuta,  Y.  Okajima and  T.  Suzuki: Scr. Mater. 55 (2006) 1095–1098. doi:10.1016/j.scriptamat.2006.08.043
• 15)   Y.  Shibuta,  Y.  Okajima and  T.  Suzuki: Sci. Technol. Adv. Mater. 8 (2007) 511–518. doi:10.1016/j.stam.2007.08.001
• 16)   A.  Roy,  H.W.  Jang and  P.-R.  Cha: Appl. Phys. Lett. 113 (2018) 193301. doi:10.1063/1.5051499
• 17)   A.  Roy and  P.-R.  Cha: J. Appl. Phys. 128 (2020) 205102. doi:10.1063/5.0026350
• 18)   A.  Roy,  M.-G.  Cho and  P.-R.  Cha: J. Appl. Phys. 131 (2022) 185106. doi:10.1063/5.0087758
• 19)   Y.  Okajima,  Y.  Shibuta,  T.  Tsuchiya,  S.  Yamaguchi and  T.  Suzuki: Appl. Phys. Express 3 (2010) 065202. doi:10.1143/APEX.3.065202
• 20)   L.  Liang,  Y.  Qi,  F.  Xue,  S.  Bhattacharya,  S.J.  Harris and  L.Q.  Chen: Phys. Rev. E 86 (2012) 051609. doi:10.1103/PhysRevE.86.051609
• 21)   L.  Liang and  L.-Q.  Chen: Appl. Phys. Lett. 105 (2014) 263903. doi:10.1063/1.4905341
• 22)   S.  Kim,  S.  Choi and  W.  Lu: ACS Nano 8 (2014) 2369–2376. doi:10.1021/nn405827t
• 23)   R.  Dong,  D.S.  Lee,  W.F.  Xiang,  S.J.  Oh,  D.J.  Seong,  S.H.  Heo,  H.J.  Choi,  M.J.  Kwon,  S.N.  Seo,  M.B.  Pyun,  M.  Hasan and  H.  Hwang: Appl. Phys. Lett. 90 (2007) 042107. doi:10.1063/1.2436720
• 24)   S.  Seo,  M.J.  Lee,  D.H.  Seo,  E.J.  Jeoung,  D.S.  Suh,  Y.S.  Joung,  I.K.  Yoo,  I.R.  Hwang,  S.H.  Kim,  I.S.  Byun,  J.S.  Kim,  J.S.  Choi and  B.H.  Park: Appl. Phys. Lett. 85 (2004) 5655–5657. doi:10.1063/1.1831560
• 25)   M.  Janousch,  G.I.  Meijer,  U.  Staub,  B.  Delley,  S.E.  Karg and  B.P.  Andreasson: Adv. Mater. 19 (2007) 2232–2235. doi:10.1002/adma.200602915
• 26)   J.J.  Yang,  M.D.  Pickett,  X.  Li,  D.A.A.  Ohlberg,  D.R.  Stewart and  R.S.  Williams: Nat. Nanotechnol. 3 (2008) 429–433. doi:10.1038/nnano.2008.160
• 27)   D.H.  Kwon,  K.M.  Kim,  J.H.  Jang,  J.M.  Jeon,  M.H.  Lee,  G.H.  Kim,  X.S.  Li,  G.S.  Park,  B.  Lee,  S.  Han,  M.  Kim and  C.S.  Hwang: Nat. Nanotechnol. 5 (2010) 148–153. doi:10.1038/nnano.2009.456
• 28)   M.J.  Lee,  C.B.  Lee,  D.  Lee,  S.R.  Lee,  M.  Chang,  J.H.  Hur,  Y.B.  Kim,  C.J.  Kim,  D.H.  Seo,  S.  Seo,  U.I.  Chung,  I.K.  Yoo and  K.  Kim: Nat. Mater. 10 (2011) 625–630. doi:10.1038/nmat3070
• 29)   G.  Bersuker,  D.C.  Gilmer,  D.  Veksler,  P.  Kirsch,  L.  Vandelli,  A.  Padovani,  L.  Larcher,  K.  McKenna,  A.  Shluger,  V.  Iglesias,  M.  Porti and  M.  Nafría: J. Appl. Phys. 110 (2011) 124518. doi:10.1063/1.3671565
• 30)   G.  Bersuker,  Y.  Jeon and  H.R.  Huff: Microelectron. Reliab. 41 (2001) 1923–1931. doi:10.1016/S0026-2714(01)00120-2
• 31)   S.  Kim,  W.  Kim and  T.  Suzuki: Phys. Rev. E 58 (1998) 3316–3323. doi:10.1103/PhysRevE.58.3316
• 32)   S.G.  Kim,  W.T.  Kim and  T.  Suzuki: Phys. Rev. E 60 (1999) 7186–7197. doi:10.1103/PhysRevE.60.7186
• 33)   S.L.  Wang,  R.F.  Sekerka,  A.A.  Wheeler,  B.T.  Murray,  S.R.  Coriell,  R.J.  Braun and  G.B.  McFadden: Phys. D 69 (1993) 189–200. doi:10.1016/0167-2789(93)90189-8