Simulation of TiN/HfO2/Pt memristor I–V curve for different conductive filament thickness

The operation of the TiN/HfO2/Pt bipolar memristor has been simulated by the finite elements method using the Maxwell steady state equations as a mathematical basis. The simulation provided knowledge of the effect of conductive filament thickness on the shape of the I–V curve. The conductive filament has been considered as the highly conductive Hf ion enriched HfOx phase (x < 2) whose structure is similar to a Magneli phase. In this work a mechanism has been developed describing the formation, growth and dissolution of the HfOx phase in bipolar mode of memristor operation which provides for oxygen vacancy flux control. The conductive filament has a cylindrical shape with the radius varying within 5–10 nm. An increase in the thickness of the conductive filament leads to an increase in the area of the hysteresis loop of the I–V curve due to an increase in the energy output during memristor operation. A model has been developed which allows quantitative calculations and hence can be used for the design of bipolar memristors and assessment of memristor heat loss during operation.


Introduction
Researchers currently show strong interest to new computer technologies such as quantum computers and neuromorphic systems. A neuromorphic system is an artificial object imitating the work of human brain. The operation principle of neuromorphic systems is "memorizing" of new information by changing the conductivity of the contacts between artificial neurons (synapses). One possible option to implement this system is a memristor array. A memristor is a device having two functional electrodes. During memristor operation its top electrode is fed with direct voltage of different sign while the bottom electrode is earthed. Transition metal oxides e.g. TiO 2 , HfO 2 , NiO, Ta 2 O 5 are typically used for the fabrication of the memristor working bodies. After the voltage is switched off the memristor does not change its state thus memorizing the last resistance value. During memristor operation its operation mode switches between the high resistance state (HRS) and the low resistance state (LRS). Memristor operation mode switching is achieved due to the formation and destruction of conductive filaments in the memristor working body. These filaments are high conductivity regions in the form of either clusters of positively charged oxygen vacancies with specific charge transfer mechanisms [1] or a separate phase having a higher conductivity as compared with that of the memristor working body [2][3][4]. Examples of phases that form the conductive filaments are the titanium ion enriched Magneli phase Ti 4 O 7 [2], the structurally ordered HfO x phase (x < 2) [3] and the TaO x phase [4].
There are no accurate data on the content of oxygen ions in the two latter phases but, by analogy with the Magneli phase, these phases are enriched with Hf [3] or Ta [4] ions, respectively. Experimental studies of the electrical conductivity of the Ti 4 O 7 phase for massive specimens have shown that a number of metal-semiconductor type phase transitions occur in this material providing for metallic conductivity type in the 150-300 K range [5,6]. As regards the HfO x and TaO x phases, their conductivity has not been studied for massive specimens. However, the formation of these phases in the conductive filaments and the higher content of metallic ions in these phases in comparison with HfO 2 and Ta 2 O 5 oxides testify to their high metallic type conductivity.
The I-V curves of bipolar memristors have hysteresis loops enabling the use of these electric devices as resistive memory cells. Different I-V curve branches correspond to different memristor operation modes, i.e., LRS and HRS. Local morphological changes occurring in TiO 2 based memristors due to electroforming were studied earlier using atomic force microscopy [7]. Experimental implementation of different I-V curve shapes recorded directly in the region to which voltage is fed allowed the authors [7] to conclude that the shape of memristor I-V curve and the type of damage produced by electrical forming correlate. It is safe to assume that along with external morphological changes in the conductive filament region caused by electroforming, the shape of memristor I-V curve is also controlled by the conductive filament thickness. However, there are no specific data on the effect of conductive filament thickness on the shape of memristor I-V curve. Nevertheless, conductive filament thickness can be critical since it largely determines the heat release conditions during memristor operation [8].
The aim of this work is to study the effect of conductive filament thickness on the shape of the I-V curve for a TiN/HfO 2 /Pt bipolar memristor by numerical simulation of memristor operation mode with the finite elements method using the Maxwell steady state equations as a mathematical basis. This approach can be referred to as "first principle" simulation of I-V curves. The Hf enriched HfO x phase (x < 2) was considered to be the conductive filament. Hafnium oxide is widely used for the fabrication of bipolar memristors in which unlike titanium oxide based memristors a wider range of electrode pairs are used, i.e., Hf-TiN [3,9], Pt-TiN [10], TiN-TiN [11,12] and Ni-TaN [13], this simplifying the choice of electrodes for the simulation of bipolar memristor operation.

Memristor model
The Maxwell steady state equations have the following form: where i is the electric current density, σ is the electrical conductivity, Е is the electric field, φ is the electric potential, ε 0 is the electric constant, ε is the relative dielectric constant and ρ is the electric charge density. The Maxwell equations allow one to calculate the electric current I in the memristor as a function of the voltage U fed to the memristor top electrode for different conductive filament heights. During the simulation the memristor was considered as a capacitor consisting of two electrodes with a HfO 2 layer located between them. The memristor model was constructed in a cylindrical coordinate system (Fig. 1). The radius r f of the HfO 2 phase layer growing in the shape of a cylinder varied in the 5-10 nm range, the radius of the entire cylindrical structure being r m = 2r f . The height of the HfO 2 layer was accepted to be 5 nm, with the conductive filament height h varying in the 0-5 nm range. The heights of the top and bottom memristor electrodes were 70 nm. When h reached 5 nm the electric circuit became short, the memristor changed to LRS and the current in the circuit obeyed Ohm's law. The memristor component dimensions accepted in this model were similar to those of actual devices [3,[9][10][11]. The conductive filament thickness in the model was varied within the dimensions of structural defects produced during memristor electroforming [2][3][4]14]. Figure 2 shows the array of finite elements constructed for the above described memristor design corresponding to the set of the Maxwell steady state equations. The behavior of single elements of the array was considered to result from linear interaction between adjacent nodes due to external forces (electric field) and be described by the respective matrix equations [15]. The minimum cell dimension was 0.2 nm in the conductive filament and 4 nm in other model regions. This choice of array parameters is stipulated by the high current density and electric field gradients that may occur at large h close to the boundary values. The finite elements array shown in Fig. 2 was constructed in the COMSOL Multiphysics software.

Choice of memristor electrode material
The material of the top memristor electrode was platinum Pt and that of the bottom electrode was titanium nitride TiN. The choice of platinum for the top electrode was stipulated by its unique properties. Under different conditions platinum can either block oxygen ions [16] or allow their passage through the electrode (i.e., be transparent) [17] thus favoring the oxidation-reduction reactions in the vicinity of the platinum/transition metal oxide interface which play an important role in memristor operation [8,17] by controlling the vacancy flux rate. Platinum has an oxidation degree of +2 (with the respective ionization energy 18.56 eV) and under a positive electric potential it enters into a chemical reaction with oxygen anions of HfO 2 , this reaction having the following form in the Kröger notation [18]: where is the site oxygen anion (in accordance with Kröger's representations [18] the site cations and anions in the lattice of an ionic crystal are in a neutral state) and is the positively charged oxygen vacancy. An important property of Pt stipulating its use as the hafnium oxide based memristor top electrode is the ability to exhibit catalytic properties, i.e., decomposition, by chemisorption, of gas molecules adsorbed on its surface. Earlier [17] it was shown that chemisorption of oxygen molecules from air occurs by the following reaction: where V ad is the vacant adsorption site in platinum and O ad is the adsorbed neutral atom (adatom) of oxygen. As a result of oxygen adatom diffusion into the platinum electrode (occurring predominantly by grain boundaries) a vacant adsorption site is left on the platinum surface. On its path the oxygen adatom traps electrons from the platinum conduction band, acquires a negative charge  and then as a result of recombination with the positively charged oxygen vacancy in the surface region becomes a neutral site anion. The reaction describing this process is as follows [17]: (4) Thus the catalytic activity of Pt causes secondary oxidation of hafnium oxide, i.e., substitution of oxygen vacancies formed by reaction (2) for oxygen ions delivered from the environment. The process of secondary oxidation is illustrated in Fig. 3. When a negative voltage is applied to the top Pt electrode, the exchange rate between oxygen adatoms having two trapped electrons and positively charged oxygen vacancies increases in accordance with electrostatic laws. Titanium nitride used as the bottom electrode material had inert properties in the model used.

Formation, growth and dissolution of the conductive filament
The crystallographic structure of the HfO x phase was not considered in the cited earlier work [3]. It is however understood that HfO 2 enrichment with Hf ions of the initial lattice can occur in two possible ways: -direct substitution of part of oxygen ions for hafnium ions; -increase in the number of oxygen vacancies. In the model used it was assumed that during the formation of the HfO x phase, the enrichment of the initial HfO 2 lattice with Hf ions occurs by the second scheme. The mechanism of oxygen vacancy formation in an ionic crystal lattice is described by the following reaction: (5) This reaction results in the delivery of two free electrons to the conduction band thus increasing the metallic type conductivity in the HfO x phase forming during memristor operation. If the number of oxygen vacancies is sufficient and the vacancy subsystem has an ordered structure (this is a distinctive feature of the Magneli phase Ti n O 2n-1 [19]), the HfO x phase should have a structure similar to that of Ti n O 2n-1 type phases thus favoring the metallic type electron conductivity mechanism.
The memristor operation was controlled by a bipolar signal having a triangular shape (Fig. 4). At Section I of the bipolar signal, oxygen vacancies permanently generated in the vicinity of the Pt electrode drift toward the region of the inert TiN electrode (in this case acting as the anode) thus providing for the formation and growth of the oxygen depleted HfO x phase. At Sections II and III of the signal until the maximum (by absolute value) negative voltage is fed to the top Pt memristor electrode, the conductive filament formed at Section I (the HfO x phase) remain stable due to the inertial character of the phase formation processes. At Section II of the signal the generation of oxygen vacancies continues, resulting in the formation of excess vacancies in the conductive filament (relative to the concentration of "site" oxygen vacancies in the HfO x phase). As a result of a change of voltage polarity at the Pt electrode the excess vacancies formed at Section II of the signal are compensated at Section III of the signal because the rate of secondary oxidation in the region close to the Pt electrode increases and additional vacancy sinks are activated in the system. At Section IV the vacancy fluxes (arising due to the continuing activity of the vacancy sinks) directed from the inert TiN electrode toward the Pt electrode trap "site" oxygen vacancies in the conductive filament, resulting in the dissolution of the HfO x phase.
The suggested mechanism of memristor operation switching to the LRS mode based on the formation and growth of the conductive HfO x phase is a heterogeneous process that includes multiple stages which may occur either in sequence or simultaneously [20]. An important feature of the heterogeneous processes is that one of the process stages is typically the limiting one [20]. Bearing in mind that vacancy drift is naturally associated with HfO x phase growth rate (in the basic case h = vt where v is the vacancy drift velocity and t is the time), it is necessary to understand under which conditions the transformation of the HfO 2 phase to the HfO x one can be disregarded during problem analysis. The theory of two-stage heterogeneous processes is the best developed. The case considered can be described with two sequential stages: oxygen vacancy drift and HfO x phase formation.
We accept that the formation of the HfO x phase in the cathode region occurs by a first order chemical reaction, i.e., HfO 2 → HfO x , whose rate ω 2 obeys the equation ω 2 = kc x where k is the reaction rate constant ([k] = с -1 ) and c x is the equilibrium concentration of oxygen vacancies in the phase HfO x expressed in molar fractions. Since the formation of the HfO x phase causes consumption of oxygen vacancies, a necessary condition for the HfO 2 → HfO x reaction to start is that the inequality c 0 > c x is met, where c 0 is the concentration of oxygen vacancies far from the interphase boundary. This condition provides for the concept of a narrow boundary layer adjacent to the interphase boundary, the rate of material supply in this layer being described by the mass transport coefficient β ([β] = с -1 ). The rate of material supply ω 1 in the boundary layer is described by the expression ω 1 = β(c 0 -c x ). If the supply of vacancies to the interphase boundary is the limiting stage of the process (β << k), then for the steady state case (ω 1 = ω 2 ) ω 2 = βc 0 [20], i.e., the formation of the HfO x phase is controlled by merely external mass transport. Thus, the rate of the HfO 2 → HfO x reaction can be ignored for the steady state case when analyzing the process, and the basic expression for I-V curve simulation will be the expression of oxygen vacancy drift velocity in the electric field. It should be noted that the accepted limitations for the formation of the HfO x phase (where the stage of vacancy supply to the interphase boundary is the limiting one and the entire process occurs in a steady state regime) are commonly realized in practice and do not contradict to the heterogeneous process nature standpoints.

I-V curve construction
The general expression of charged vacancy drift velocity v in the electric field Е is as follows [21]: where a is the crystal lattice spacing, D v is the oxygen vacancy diffusion coefficient, q is the oxygen vacancy charge, k B is the Boltzmann constant and T is the Kelvin temperature. Introducing the characteristic electric field magnitude E 0 = 2k B T/qa and using for a vacancy drifting in the electric field such a kinetic constant as the mobility m v (m v = qD v /k B T), one can rewrite Eq. (6) in the following form (7) which is suitable for simulation. For the migration of oxygen vacancies in hafnium oxide at room temperature (T = 300 K) E 0 = 5 ⋅ 10 7 V/m which for the memristor model used corresponds to the voltage U 0 = 0.25 V. E 0 is typically used for defining the electric fields at which an electric device operates. The magnitude of weak electric fields E is far lower than E 0 (E << E 0 ) whereas fields for which the ratio E ≈ E 0 is met should be considered strong [22]. For the simulation used, Eq. (7) should be rewritten to correlate v with the variable U.
The expression of the HfO x phase layer height as a function of time (h = vt) is true for constant voltage U. If h is represented as the two-parameter function h = h(U, t) (which is the case in this simulation), taking into account that y = sinh (U/U 0 ) is an exponential-type function, the expression h = h(U, t) can be represented as follows: (8) where K 0 is a semi-empirical constant which takes into account the typical duration of the memristor operation cycle. The use of Eq. (8) for I-V curve calculation allows one to exclude the time parameter from problem consideration.
The typical voltage of switching to LRS mode for bipolar memristors is 1.0-1.5 V [3,[9][10][11][12]. In this simulation the switching voltage is limited to 1.0 V. Substitution of h = 5 nm and U = 1.0 V into Eq. (8) yields K 0 = 0.183 nm. This K 0 was used for the calculation of the current I at all U and for any linear memristor dimensions (with the accepted limitations). Data on the electric and dielectric properties of the materials (TiN, Pt, HfO 2 и HfO x ) which determine the design of the memristor [23,24] are shown in the Table. The simulated I-V curves for the conductive filaments of different diameters are shown in Fig. 5. It can be seen from Fig. 5 that the I-V curves have hysteresis loops. At Sections II and III of the signal, the presence of metallic conductivity type conductive filaments in the memristor structure provided for adherence to Ohm's law. The I-V curve evolution pattern at Sections I and IV of the signal is close to exponential. For r f = 5 nm the I-V curve has the maximum matching current of 1.5 mA and for r f = 10 nm this value is 45 mA. At different conductive filament thicknesses, application of the same signal generates different currents in the memristor, both in the LRS and in the HRS regimes. Thus the simulated I-V curve is sensitive to the conductive filament thickness and hence the model developed in this work allows one to take into account not only the direct correlation between the current and the voltage but also the effect of the conductive filament thickness on the current in the memristor. The latter fact is important for the design of memristors bearing in mind the heat loss during memristor operation. HfO x 2 · 10 4 * -10 6 Pt 5 · 10 6 -10 6 * σ for the HfO x phase was calculated based on earlier I-V curve [9] and in accordance with earlier conductive filament dimensions data [3].

Conclusion
The operation of the TiN/HfO 2 /Pt bipolar memristor was simulated by the finite elements method using the Maxwell steady state equations as a mathematical basis. The conductive filament was the HfO x phase possessing metallic conductivity type. The memristor operation mode included four sequential time intervals corresponding to different sections of the bipolar signal having a triangular shape. Simulation was carried out for the conductive filament of a cylindrical shape whose radius was varying within 5-10 nm. The simulated I-V curves of the memristor had hysteresis loops which fact agrees with earlier experimental data. Greater conductive filament thickness corresponded to wider hysteresis loops testifying to a higher energy output during memristor operation. At sections I and IV of the bipolar signal the I(U) dependence was exponential while at Sections II and III Ohm's law was obeyed. Depending on conductive filament thickness the memristor current varied in the same voltage range not only in LRS mode but also in HRS one. The I-V curve calculation method described in this work can be used for analyzing heat loss during memristor operation.