Non-linear memristor switching model

We introduce a thermodynamical model of filament growing when a current pulse via memristor flows. The model is the boundary value problem, which includes nonstationary heat conduction equation with non-linear Joule heat source, Poisson equation, and Shockley- Read-Hall equations taking into account strong electron-phonon interactions in trap ionization and charge transport processes. The charge current, which defines the heating in the model, depends on the rate of the oxygen vacancy generation. The latter depends on the local temperature. The solution of the introduced problem allows one to describe the kinetics of the switch process and the final filament morphology.


Introduction
In recent years significant efforts are being made for implementation of fourth electronic elementresistor with memory called memristor. This element is a passive device whose resistance can be changed by passing through a certain electric current pulse. The memristor has a simple metal-high-kmetal structure.
Functionally the memristors can be divided into two-level (one-bit) and multilevel (multibit) ones. If the one-bit device may be used as a cell of non-volatile memory with a large storage time and high resistance to radiation, the multibit memristors on the one hand can increase the density of integration of the non-volatile memory, and on the other hand can be used to form the basis of adaptive neuromorphic (cognitive) computational systems. One of the most interesting and hitherto unsolved problems is the nature of the intermediate (discrete) multi-level states in a memristor.
The hypothesis that the effect of memristive switch is directly connected with the formation of a filament, which occurs when certain current pulse passes over a dielectric, is distributed among specialists. However, the kinetics of filament formation during the transition from the high resistance state of the memristor to a low resistive state is still not clear. Another unsolved problem in the way of developing memristor arrays is a description of the forming process of memristor elements -the first switch from the initial high resistance state to a low resistance at high voltages, accompanied by intense Joule heat release. Nature of the forming of memristor matrix is still unclear, and the absence of any plausible physical model constrains the practical implementation of these elements.
It is considered that charge transport in high-κ is associated with the presence of defects. The most probable and the most widespread defects are oxygen vacancies. Properties of these vacancies are actively studied both experimentally and theoretically. However, many contradictory data can be found in literature, and the role and characteristics of oxygen vacancies in the charge transport through dielectric process are not established yet. It is obvious that rate of change in structure of dielectric (namely the generation rate of oxygen vacancies associated with conductivity) in strong electric field is mainly determined by non-stationary processes of heat and mass transfer. Thus, the kinetics of such processes determines final morphology and properties of memristor.

The theory
Physics model of resistive switch in thin films of oxides of transition metals consists of 4 main subtasks: strongly nonlinear and non-stationary thermophysical problem, problem of finding charge distribution and charge transfer, problem of finding of the heat source function and problem of generation, recombination and migration of defects that act as traps for electrons and holes.

Thermophysical problem
Generation (recombination) of defects in active medium of ReRAM element is assumed to be induced by local heat release caused by electric current flowing through it. For modeling purposes the structure of memristor is assumed to consist of 3 layers: 2 half-infinite layers of metal electrodes with dielectric layer between them as shown in Fig 1. Layer thickness is d. Temperature field in the volume of memristor is described by heat conduction equation which we regard as having cylindrical symmetry: (1) where ρ is the medium density, c is the medium heat capacity, T is the medium local temperature, t is time, r and z are radius and height, respectively; λ is heat conductivity, Q is heat per unit volume (will be defined later on).
The boundary conditions are as follows: no heat flux on the axis of symmetry; the temperature far away from filament equals the initial temperature T 0 ; heat fluxes and temperatures at the dielectricelectrode border are equal.

Charge distribution and charge transfer problem
In solid bodies the charge carriers are electrons and electron holes -quasi-particles that represent collective reaction of electrons on valence atomic orbitals of elements that form the medium [1,2]. The distinctive feature of such quasi-particles is that they have a charge equal to the charge of free electron q in absolute value and that they have an effective mass m * defined by crystal structure of the substance. It is generally accepted that the charge of electron hole is positive and the charge of electron is negative. Later on we will assume that elementary charge q is positive.
Charge transport involves both «free» electrons and holes (1 in conduction band and 1' in valence band), and trapped quasi-particles (Fig. 2). The trapped quasi-particles can be ionized in respective zone i.e. become «free» (2, 2' in Fig. 2), «free» charge carriers can become trapped (3, 3' in Fig. 2), providing the closeness of traps trapped electrons and holes can tunnel between adjacent traps (4, 4' in Fig. 2). Also free electrons can recombine with trapped holes (5 in Fig. 2), free holes can recombine with trapped electrons (5' in Fig.2). Recombination of free electrons and free holes can be neglected because of low concentration of quasi-particles in carrier band and valence band of the dielectric. Abovementioned processes are described by Shockley-Read-Hall equations complemented by effect of strong electron-phonon interaction in the processes of ionization of traps and charge transport: where n is concentration of free electrons, v n is drift velocity of free electrons, e = F/F is a unit vector, F is electric field, σ n is capture cross-section of electron on the trap, N is concentration of charge traps (assuming that the trap for electron and a trap for electron hole are the same defects, namely, oxygen vacancies), n t is concentration of trapped electrons, P ion is ionization rate of trapped charge carriers, σ r is recombination cross-section for electron-hole pair, p t is concentration of trapped holes, p is concentration of free holes, v p is drift velocity of free holes, σ p is capture cross-section of hole on trap. Terms (∂n t /∂t) tun and (∂p t /∂t) tun represent charge transfer caused by tunneling between traps: (6) where s = N -1/3 -average distance between traps, P tun is rate of tunneling between traps.
Process of ionization of trapped into conduction/valence band electrons/holes is described by the theory of multiphonon ionization of neutral isolated trap (2, 2' in Fig. 2): where W ph is phonon energy (breathing mode), W t is thermal trap energy, W opt is optical trap energy, k is the Boltzmann constant, I n are modified Bessel functions, ћ = h/2π is the Plank constant, P i is probability of tunneling through the barrier W [3].
Jumping between adjacent traps is described by phonon-assisted tunneling between traps: where [4]. Filling of the first trap derives from contact and described by the following equations: (10) where coefficients of filling of the trap with electrons (holes) from the contact v inj and back ionization of the trap into contact v ion were found using results of Refs. [4,5]. Ionization rate for phonon-bounded trap is described by following equation: Here Φ is energy difference between the energy of the bottom of conductive band (for electrons, the top of the valence bands for holes) of the dielectric and Fermi level of the contact, a is the distance from the trap to the contact, V out is velocity of released from the trap electron (hole) in the contact. Injection velocity v inj from contact to the trap can be found using the fact that in stationary conditions trap's filling complies statistics with Fermi level in the contact: Relationship between injection velocity and ionization velocity is found from stationary equations (9), (10) without charge transfer to adjacent traps: (13) Filling of the last trap can be calculated similarly with filling of the first with substitution of F with -F. Also, it is necessary to qualify that for electrons the first trap is the trap near the electrode with higher potential, for the holes -with lower one.
Total current j includes electron j n and hole j p components: where μ ~ (0,1-10) cm 2 /V·s is charge carrier mobility in dielectric layer.
It is necessary to note that free and trapped charge carriers in the volume of dielectric layer distort the external electric field, and their influence needs accounting for by solving Poisson equation: where φ is electric field potential, ε is low-frequency permittivity of the medium.

Heat source function
Heat source function is caused by the release of heat due to inelastic scattering of free charge carriers (Joule heat) and recombination of electrons and holes. Generally speaking, the process of recombination can go with creation of photons or phonons. In the process of recombination of free electron with trapped hole (or free hole with trapped electron) the full momentum is nonzero and remains so after the recombination. It is known that photons have negligible momentum and phonons can have momentum comparable to ones of electrons and holes. So, in this model creation of photons can be considered highly improbable process and ignored. Thus, in recombination of opposite charges all energy releases in the form of phonons and should appear in the heat source function: Here E g is the band-gap energy of dielectric, W t n and W t p are thermal energies of electron and hole traps, respectively. The first member in Eq (19) represents Joule heat, the second one describes recombination of free electrons on trapped holes, the third one is recombination of free holes on trapped electrons. It is noteworthy that for values characteristic for this problem: electric field ~10 6 V/cm and current 10 -6 A/cm 2 , Joule heat have the value ~1 W/cm 3 . For values σ r ~ 10 -13 cm 2 , v n,p ~ 10 7 cm/s, n, p ~ 10 7 cm -3 , n t , p t ~ 10 18 cm -3 , E g ~ 6 eV, W t ~ 1 eV released energy from recombination can be ~(0.1÷10) W/cm 3 . So, the heat release in recombination in this case cannot be neglected! The change in energy of charge carriers in tunneling between traps is taken into account in Eq (19) using equations (14)-(16).

Generation, recombination and migration of oxygen vacancies
As was mentioned above forming and resistive switch to low-resistance state is related to trap generation for charge carriers. Oxygen vacancies are acting as traps. In accordance with experimental data certain fraction of vacancies exists in dielectric immediately after synthesis. However, defects generated during forming are responsible for performance attributes. Switching to high-resistance state is caused by recombination and spatial redistribution of defects described by following system of equations: where N 0 is atomic concentration, s 0 = N 0 -1/3 is interatomic distance, W ox is energy of interstitialvacancy pair formation, N ox is concentration of interstitial atoms of oxygen created in formation of interstitial-vacancy pair, v ox is drift velocity of interstitial atoms of oxygen, σ ox -recombination crosssection for interstitial-vacancy pair. The first term of the right side of (20) describes the generation rate for new interstitial-vacancy pair with lowering of energy of formation of the pair in strong electric field [6]. As a pre-exponential factor the concentration of oxygen atoms in their atom positions is used, frequency factor W ph /h is calculated using characteristic oscillation frequencies of atoms near their equilibrium position i.e. phonon oscillation modes. The second term in (20) describes the recombination of oxygen vacancies with interstitial oxygen. The third one is associated with the drift of oxygen vacancies occurring as jumping of adjacent oxygen atoms to the vacancy position by Hill mechanism [7], at that initial vacancy recombines and in a place of atom new vacancy forms. Equation (21) describes generation/recombination of interstitial atoms of oxygen. The first and the second terms in (21) are similar by design with respective members of Eq. (20). The third one describes the drift of interstitial oxygen in electric and thermal fields. Dependence of drift velocity v ox on temperature of the medium and electric field requires additional investigation.
To find the recombination cross-section of interstitial-vacancy pair following considerations were used. If interstitial oxygen atom gets to the area near the vacancy with attractive potential the interstitial-vacancy pair recombines. Within Frenkel model cross-section of such area is equal to the area of a circle with radius matching the distance between the center of the vacancy and the maximum of the potential: Equation (22) shows that area of recombination cross-section decrease with increase of electric field (voltage on the structure) for interstitial-vacancy pair. Thus, for linear dependence of drift velocity of interstitial-vacancy pair from electric field the probability of recombination does not depend on electric field. However, probability of generation grows exponentially with increase of electric field and temperature. So, in current problem the switching from high-resistance state to low-resistance state is anticipated in the area of local heating. At that the electric field will promote the transport of formed interstitial sites from the area of filament and migration of vacancies of oxygen to the active area. Under the electric field of opposite direction, the interstitial sites will migrate to the active area of filament and recombine with vacancies of oxygen. In the process oxygen vacancies will leave the area of filament. In conjunction the processes will lead to the switch to high-resistance state.

Summary
In conclusion, this study introduces thermodynamical model of filament growing when a current pulse of variable length and value flows. The model is the boundary value problem, which includes nonstationary heat conduction equation with non-linear Joule heat source, Poisson equation, and Shockley-Read-Hall equations taking into account strong electron-phonon interactions in trap ionization and charge transport processes. The charge current, which defines the heating in the model, depends on the rate of the oxygen vacancy generation. The latter depends on the local temperature. Kinetic parameters of our model are obtained from the transport and optical experiments and ab initio simulations. The solution of the introduced problem allows us to describe the kinetics of the switch process and the final filament morphology. The resulting filamentary structure in turn may explain the memristor switch.
End-to-end solution of the task allows us to develop recommendations for optimizing the memristor technology and predict the behavior of memristor arrays, including neuromorphic electronics and non-volatile large memory arrays