An Enhanced Lumped Element Electrical Model of a Double Barrier Memristive Device

TL;DR

This paper presents an improved lumped element electrical model for a double barrier memristive device, enabling fast and physically meaningful simulations suitable for neuromorphic circuit applications.

Contribution

The paper introduces a novel enhanced electrical model for a VCM-based double barrier memristive device that balances simulation speed and physical accuracy.

Findings

Model enables rapid circuit simulations

Verification against kinetic Monte-Carlo simulations

Maintains physically meaningful parameters

Abstract

The massive parallel approach of neuromorphic circuits leads to effective methods for solving complex problems. It has turned out that resistive switching devices with a continuous resistance range are potential candidates for such applications. These devices are memristive systems - nonlinear resistors with memory. They are fabricated in nanotechnology and hence parameter spread during fabrication may aggravate reproducible analyses. This issue makes simulation models of memristive devices worthwhile. Kinetic Monte-Carlo simulations based on a distributed model of the device can be used to understand the underlying physical and chemical phenomena. However, such simulations are very time-consuming and neither convenient for investigations of whole circuits nor for real-time applications, e.g. emulation purposes. Instead, a concentrated model of the device can be used for both fast simulations and real-time applications, respectively. We introduce an enhanced electrical model of a valence change mechanism (VCM) based double barrier memristive device (DBMD) with a continuous resistance range. This device consists of an ultra-thin memristive layer sandwiched between a tunnel barrier and a Schottky-contact. The introduced model leads to very fast simulations by using usual circuit simulation tools while maintaining physically meaningful parameters. Kinetic Monte-Carlo simulations based on a distributed model and experimental data have been utilized as references to verify the concentrated model.

Figures

Equations

$$k=\nu {\text{rme}}^{-\frac{{\Phi }_{\mathrm{a}}}{k_{\mathrm{B}}\vartheta }}.$$

$$u(t)$$

$$\dot{\mathbf{z}}$$

$$z=\frac{x-x_{\min }}{x_{\max }-x_{\min }},\mathrm{with}z\in [0,1],$$

$$\dot{z}=\frac{-\dot{Z}w(z)}{e^{{\phi }_{\mathrm{a}}(u,z)}}\sinh \left(\frac{u_{\mathrm{r}}(u,u_{\mathrm{s}},z)+u_{\mathrm{e}}-U_{\mathrm{c}}}{U_{\mathrm{e}}}\right)$$

$$\dot{Z}=2\nu a,\mathrm{with}a=\frac{d_{\mathrm{hop}}}{x_{\max }-x_{\min }}$$

$$U_{\mathrm{e}}=\frac{2}{c}\frac{d_{\mathrm{e}}}{d_{\mathrm{hop}}}{U}_{\vartheta },\mathrm{with}{U}_{\vartheta }=\frac{k_{\mathrm{B}}\vartheta }{q_{\mathrm{e}}}.$$

$$u_{\mathrm{r}}(u,u_{\mathrm{s}},z)=\sigma (-u)\left[1-z\right]u_{\mathrm{s}},$$

$$\sigma(\xi)=\left\{\begin{array}[]{rcl}1&\mathrm{for}&\xi>0\\ 0&\mathrm{otherwise}&\end{array}\right.$$

$$w(z)=\left[1-2w_0\right]\left[1-{\left[2z-1\right]}^{2p}\right]+w_0,$$

$${\phi }_{\mathrm{a}}(u,z)=\sigma (u)\left[{\phi }_{{\mathrm{a}}_1}+z\left[{\phi }_{{\mathrm{a}}_0}-{\phi }_{{\mathrm{a}}_1}\right]-{\phi }_{{\mathrm{a}}_{\mathrm{r}}}\right]+{\phi }_{{\mathrm{a}}_{\mathrm{r}}},$$

$$i_{\mathrm{s}}(u_{\mathrm{s}},z)=I_{\mathrm{s}}{\text{rme}}^{-\left[{\phi }_{\mathrm{s}}(z)+{\alpha }_{\mathrm{f}}\sqrt{\frac{|u_{\mathrm{s}}|-u_{\mathrm{s}}}{{\alpha }_{\mathrm{s}}{U}_{\vartheta }}}\right]}\left[{\text{rme}}^{\frac{1}{n(z)}\frac{u_{\mathrm{s}}}{U_{\vartheta }}}-1\right].$$

$${\alpha }_{\mathrm{s}}=\frac{2d_{\mathrm{s}}}{D_{\mathrm{s}}},\mathrm{with}{D}_{\mathrm{s}}=\frac{q_{\mathrm{e}}^2}{4\pi {ϵ}_0{ϵ}_{\mathrm{r}}{k}_{\mathrm{B}}\vartheta },$$

$${\phi }_{\mathrm{s}}(z)={\phi }_{{\mathrm{s}}_0}+z\left[{\phi }_{{\mathrm{s}}_1}-{\phi }_{{\mathrm{s}}_0}\right]\mathrm{and}$$

$$n(z)=n_0+z\left[n_1-n_0\right],$$

$$u_{\mathrm{s}}={\hat{R}}_{\mathrm{s}}(u_{\mathrm{s}},z)i_{\mathrm{s}}(u_{\mathrm{s}},z),$$

$${\hat{R}}_{\mathrm{e}}(z)=R_{{\mathrm{e}}_0}+z\left[R_{{\mathrm{e}}_1}-R_{{\mathrm{e}}_0}\right],$$

$$u_{\mathrm{e}}={\hat{R}}_{\mathrm{e}}(z)i_{R_{\mathrm{e}}}.$$

$$\begin{gathered} \text{eqalign}{i}_{\mathrm{t}}(u_{\mathrm{t}},z)=I_{\mathrm{t}}\frac{g(-u_{\mathrm{t}},z)-g(u_{\mathrm{t}},z)}{{\alpha }_{\mathrm{t}}^2(z)},\mathrm{with} \\ g(u_{\mathrm{t}},z)={\phi }_{\mathrm{t}}(u_{\mathrm{t}}){\text{rme}}^{-{\alpha }_{\mathrm{t}}(z)\sqrt{{\phi }_{\mathrm{t}}(u_{\mathrm{t}})}},\mathrm{where} \\ {\phi }_{\mathrm{t}}={\phi }_{{\mathrm{t}}_0}+\frac{1}{2}\frac{u_{\mathrm{t}}}{U_{\vartheta }}\mathrm{and}{I}_{\mathrm{t}}=\frac{A}{D_{\mathrm{t}}^2}\frac{k_{\mathrm{B}}{q}_{\mathrm{e}}}{2\pi \mathrm{h}}\vartheta . \end{gathered}$$

$${\alpha }_{\mathrm{t}}(z)={\alpha }_{{\mathrm{t}}_0}+z\left[{\alpha }_{{\mathrm{t}}_1}-{\alpha }_{{\mathrm{t}}_0}\right],\mathrm{with}{\alpha }_{{\mathrm{t}}_0}\le \alpha (z)\le {\alpha }_{{\mathrm{t}}_1}.$$

$$u_{\mathrm{t}}={\hat{R}}_{\mathrm{t}}(u_{\mathrm{t}},z)i_{\mathrm{t}}(u_{\mathrm{t}},z),$$

Full text

Kinetic Monte-Carlo simulations based on a distributed model of the device can be used to understand the underlying physical and chemical phenomena. However, such simulations are very time-consuming and neither convenient for investigations of whole circuits nor for real-time applications, e.g. emulation purposes. Instead, a concentrated model of the device can be used for both fast simulations and real-time applications, respectively. We introduce an enhanced electrical model of a valence change mechanism (VCM) based double barrier memristive device (DBMD) with a continuous resistance range. This device consists of an ultra-thin memristive layer sandwiched between a tunnel barrier and a Schottky-contact. The introduced model leads to very fast simulations by using usual circuit simulation tools while maintaining physically meaningful parameters.

Kinetic Monte-Carlo simulations based on a distributed model and experimental data have been utilized as references to verify the concentrated model.

• January 2017

Keywords: memristive devices, resistive switching, neuromorphic circuits, electrical modeling, nanoelectronics

1 Introduction

Resistive switching devices are essential components for today’s nonvolatile memory applications. In general, they are built by a capacitor-like metal-insulator-metal structure. Depending on material compositions, different chemical and physical effects lead to a change of the total resistance. It should be stressed that all resistive switching devices can also be interpreted as memristive systems [1], which are in general nonlinear resistors with memory [2, 3]. Especially memristive systems with a continuous resistance range are potential candidates for neuromorphic circuits [4, 5].

For our approach, we focus on a double barrier memristive device (DBMD) [6], which consist of an ultra-thin memristive layer sandwiched between a tunnel barrier and a Schottky-like contact. Several benefits as a continuous resistance range, an intrinsic current compliance, an improved retention, no need for an electric forming procedure and low power consumption make the DBMD particularly suitable for neuromorphic circuits. Just as for other electronic components, a parameter spread especially for memristive devices fabricated in nanotechnology is unavoidable. This parameter spread aggravates reproducible analyses. A kinetic Monte-Carlo model of this device with distributed parameters can help to understand the underlying chemical and physical phenomena [7, 8]. However, investigations with simulations based on a distributed model are very time-consuming. Because of this, a distributed model is neither convenient for simulations of whole circuits nor for real-time applications, e.g. emulation purposes.

Our intention is to build a replica of the device using an enhanced lumped element electrical (concentrated) model for fast simulations of both, the device itself and complex circuits including such devices. This approach also allows for reproducible analyses [9].

Based on a distributed model from [8] and on the electrical model from [6] an enhanced electrical model of the DBMD with concentrated parameters has been built up and simulated using LTSpice [10], cf. figure 1.

We have verified the concentrated model by comparisons with kinetic Monte-Carlo simulations as well as measurements. Additionally, Ochs et al[11, 12] has shown that the concentrated model can be used for emulation purposes in real-time applications.

2 The Double Barrier Memristive Device

The DBMD has been introduced by Hansen et al[6]. It consists of an ultra-thin niobium oxide NbxOy followed by an aluminum oxide Al2O3 layer sequence sandwiched between a gold Au and an aluminum Al electrode, see figure 2.

A metal-semiconductor junction at the Au/NbxOy interface creates a Schottky-contact, whereas the Al2O3 layer is an electrically high quality tunnel barrier. For the DBMD the resistance change originates from oxygen diffusion caused by an externally applied electrical field resulting in modifications of local states within the NbxOy solid state electrolyte [6, 8]. These modifications influence the interface properties of both the Schottky-contact and tunnel barrier, simultaneously, which in turn leads to an overall resistance change.

Measured Hysteresis Curve

The DBMD has been characterized by applying a triangle-shaped voltage according to the inset of figure 3 and measuring the current, normalized for a cross sectional area of $1~{}\mu\mathrm{m}^{2}$, cf. [6] and [8]. When plotting the current versus the voltage, see figure 3, the familiar hysteresis curve appears that is typical for memristive devices.

A logarithmic scale increases the visibility of a continuous transition between high and low resistance state. The experimentally measured data is utilized as a reference for verifying both the distributed as well as the concentrated model.

3 Description of Device Physics

A distributed model of the device is useful for a correct interpretation of physical and chemical phenomena within the device. Dirkmann et al[8] has used such a model to investigate physical and chemical effects, which are responsible for the memristive functionality. There, investigations have been done by using a kinetic Monte-Carlo model. Based on investigations of [8] and actual insights, we have extended the known model of Hansen et al[6] to the presented concentrated model. For this progress, a qualitative recapitulation of the resistive switching behavior on the atomic level is presented in the following. A long time scale investigation with the step response of the device belongs to novel approaches within the frame of this work.

Resistive Switching Behavior

Investigations of [8] yield that the current through the device depends on the defect distribution within the memristive NbxOy layer. More precisely, the defect distribution is assumed to influence the effective thickness of the tunnel barrier as well as the Schottky-barrier height. Regarding the sign of an externally applied voltage, different physical processes occur. The following description is segmented into regimes based on the excitation of figure 3 (inset).

Thermodynamical equilibrium:

Without an external voltage, the device is in its thermodynamical equilibrium. Due to the Coulomb potential, positive and negative charges are uniformly distributed within the electrolyte, see figure 4.

This is also the high resistance state of the memristive device, because of a large Schottky-barrier height and tunnel barrier thickness.

Positive applied external voltage $0<t<T/2$:

A positive applied voltage leads to a voltage drop over the NbxOy layer. The resulting electrical field exerts a motion of charged defects. Physically, the ion motion is based on an ion-hopping phenomenon [13]. The rate of ion motion from one stable position to the other can be described by the Arrhenius-law

[TABLE]

Here $\Phi_{\mathrm{a}}$ is the activation energy of defect motion within the electrolyte in J, $\vartheta$ is the temperature in K, $\nu$ is the hopping frequency in Hz and $k_{\mathrm{B}}$ is the Boltzmann constant in J/K.

The ion motion due to the applied electrical field results in an enrichment of negative ions at the Schottky-interface and positive ions at the tunnel barrier. Regarding the interface potential at the Schottky-contact, the barrier height decreases with increasing number of ions at the contact. On the other hand, the local vacancy concentration at the Al2O3 interface is of particular importance for the resulting electron tunneling current. An increase of the vacancy concentration decreases the effective tunnel barrier thickness and therefore the tunneling current increases, cf. figure 5. In total, the device changes gradually from a high to a low resistance state.

Negative applied external voltage $T/2<t<T$:

For negative applied voltages, the Schottky-contact almost totally blocks the current, so that the voltage drops nearly completely across this contact. This results in an increased electric field near the interface. If the field becomes high enough, the oxygen ions are detached (desorbed) from the interface and move back into the electrolyte to finally reconstitute thermodynamical equilibrium. This phenomenon, which is the reset process, is illustrated in figure 6.

The desorption mechanism of ions at the Schottky-interface has to be taken into account in the concentrated model, because it dominates the time scale for negative applied voltages [8].

Long Time Scale Investigation

Within the frame of this work, firstly a step response of the DBMD device was investigated. Thereby it turned out, that the assumption of a fixed number of defects within the NbxOy electrolyte might not be valid for long time scales. In figure 7 the simulated currents assuming a constant defect number and defect formation within the electrolyte are compared with the measured current. The simulation with a fixed number of defects was done using the model presented by Dirkmann et al[8]. A constant voltage of $2.5\>\mathrm{V}$ has been chosen for the simulation as well as for the experiment. When a constant set voltage, which is sufficient to move point defects within the electrolyte, is applied to the device, the defects move towards the Au electrode. During this process, the device resistance decreases and therefore the current through the device increases, see figure 7. For longer times, all defects concentrate at the Au electrode. A further resistance change is inhibited and the current goes into saturation. In contrast to this, the measured current rises further for long time scales.

It is well-known, that a strong electric field applied at an oxide material leads to formation of point defects, as oxygen vacancies and oxygen interstitials [14]. Thus it stands to reason, that for long time scales a formation of point defects within the electrolyte region might also occur within the DBMD. This assumption has been tested by an extended version of the model presented in [8] towards the defect formation. The formation of point defects is included using the rate equation (1), see [15, 16]. The simulation result changes then with a variable number of defects within the electrolyte, cf. figure 7. The good agreement between simulation and experiment indicates a defect formation process within the DBMD for long time scales and high electric fields. It is notable, that this process is a long time process that does not affect processes on a short time scale. Although this model is able to explain the measured current for long time scales very well, it needs to be said, that it cannot be excluded that another explanation model for the current behaviour on long time scales can be found. The following kinetic Monte-Carlo simulations use a variable defect number within the electrolyte in order to verify the lumped element model.

4 Electrical Description with Concentrated Parameters

The distributed model leads to a better understanding of underlying physical and chemical effects. But for circuit simulations a concentrated model is more appropriate. In order to preserve the relationship to physical phenomena, we derive an enhanced electrical model with physically meaningful concentrated parameters. In the following, a derivation of the model and a classification as a memristive system are presented.

Memristive Systems

The resistive switching device is in general a memristive system with an internal state. A mathematical classification of the device into the theory of general memristive systems as given by Chua et al[3] is missing. The concentrated model enables integrating the DBMD into the general memristive system theory. For this reason, we start with a brief review of these systems.

The $n$th-order voltage-controlled memristive system

[TABLE]

is made up of an algebraic equation for the input-output relation together with a differential equation, which describes the memristive behavior or dynamics of the system [3]. Actually, the memristive system $\hat{R}$ is a generalized response, which interrelates the input voltage $u$ with the output current $i$. The current and voltage are scalars and hence $\hat{R}$ is a scalar function, which is generally nonlinear. A continuous nonlinear vector function $\boldsymbol{f}$ describes the dynamics of the state variable $\boldsymbol{z}$, with dimension $n$.

Concentrated Model of the DBMD

The topology of the circuit proposed by Hansen et al[6] has been used as a starting point for the electrical representation of the concentrated model. There, a model based on experimental results was developed initially. In contrast to [6] we have modified the components in the equivalent circuit with enhanced functionalities representing physical properties.

For the concentrated model, a deeper modeling approach depending on already known as well as novel physical insights is desired. On the other hand, the underlying physical and chemical phenomena are complex. Concerning the concentrated model, a preferably straightforward implementation of the functionality, without sacrificing the physical interpretation, is beneficial. Furthermore, it is advised to reformulate the mathematical description of the device in order to classify it as a memristive system. To this end, we focus the perspective rather on an electrical point of view than on a physical one by using normalized parameters and ordinary differential equations. This leads to the concentrated model depicted in figure 8, where the Schottky-contact is modeled by a diode, and the electrolyte region and the tunnel barrier are represented each by a parallel connection of a resistor and a parasitic capacitor.

However, the resulting parameters and equations, which are electrically interpretable, still have a corresponding physical meaning. In the following, physically motivated parameters and equations are deduced.

Memristive Behavior

As a first step, we derive a concentrated formulation for the memristive behavior including the internal state, cf. equation (2b). Physically, this internal state is represented by the distribution of ions in the electrolyte, which determines a total resistance value [6, 8]. Individual regions as Schottky-barrier, solid state electrolyte and tunnel barrier are coupled through this internal state. In the distributed model the position of each single ion affects the overall resistance inherently. Instead, in the concentrated model, the ion distribution is expressed by the average ion position, see figure 8, or more precisely, by the center of gravity of oxygen ions $x$. The internal state of the concentrated model is thus defined by a normalized average ion position

[TABLE]

where $x_{\mathrm{min}}$ and $x_{\mathrm{max}}$ are the minimal and maximal absolute average positions, see figure 8 (blue, dashed line). Thus, $z=0$ corresponds to the low Ohmic state, whereas $z=1$ is related to the thermodynamical equilibrium or high Ohmic state. The thermodynamical equilibrium is given for homogeneously distributed positive and negative ions. This state is represented by an average ion position right in the middle of the electrolyte. With this, we set $x_{\mathrm{min}}=0$ and $x_{\mathrm{max}}=d_{\mathrm{e}}/2$, where $d_{\mathrm{e}}$ is the electrolyte width.

Inspired from the results of Meyer et al[13], we describe the motion of ions by an ion hopping phenomenon in terms of a jump over a potential barrier $\Phi_{\mathrm{a}}$, namely activation energy. It is a material property, which together with equation (1) leads to an overall ion drift velocity. This drift velocity depends on the applied electrical field. In the concentrated model the electrical field is transformed to a corresponding voltage. Due to this transformation, nonlinearities between electrical fields and resulting voltages are neglected. With this, we can describe the normalized drift velocity of ions

[TABLE]

in terms of voltages instead of electrical fields. This velocity in turn is an ordinary differential equation describing the memristive behavior, cf. equation (2b). The drift amplitude

[TABLE]

contains the hopping frequency $\nu$ and a normalized hopping distance $a$, where $d_{\mathrm{hop}}$ denotes the absolute hopping distance per jump.

A sum of voltages causes the ion motion, see equation (2d) and figure 8. Contributing voltages are the electrolyte voltage $u_{\mathrm{e}}$, a state-dependent amount of the voltage drop over the Schottky-barrier $u_{\mathrm{r}}(u,u_{\mathrm{s}},z)$, with Schottky-barrier voltage $u_{\mathrm{s}}$ and a voltage $U_{\mathrm{C}}$ representing the Coulomb potential. The reference electrolyte voltage is defined as follows:

[TABLE]

Here, $c$ is the charge number of the used material (NbxOy), $U_{\vartheta}$ is the thermal voltage and $q_{\mathrm{e}}$ is the elementary charge.

For the reset process, a particular contribution to the sum of voltages affecting the ion motion is given by an amount of the Schottky-voltage during the reset process. This amount increases especially for ions allocated next to this interface. Therefore a state as well as voltage-dependent function modeling this amount of contribution is needed. To this end, we define a function

[TABLE]

where the dependency with respect to the state variable is assumed to be linear. Ions located directly at the interface, $z=0$, are affected by the total amount of the Schottky-barrier voltage, whereas for $z=1$ the amount of the Schottky-barrier voltage influencing the ion movement is zero. In equation (2g), the unit step function

[TABLE]

ensures that this amount of voltage applies only during the reset process, i.e. for $u<0$.

It is common practice to restrict the normalized state variable between 0 and 1 by using a window function, as shown by Biolek et al[17]. For our approach, we need a modified window function, which overcomes the boundary lock problem. Preparing the window function of [17] leads to an enhancement

[TABLE]

where the offset $w_{0}$ is positive and can be chosen arbitrarily small. This window function still involves the boundary caused nonlinearities, where $p$ controls the edge steepness.

Adsorption and desorption of ions at the Schottky-interface are further results from investigations with the distributed model [8]. Besides, from long term investigations within the scope of this work, we assume that the number of defects is not constant. For a consistent electrical model, these effects must be considered. Due to the charge separation during the set process, the activation energy increases for $z\rightarrow 0$. This is because more and more ions are impounded at this interface for positive voltages. In addition, also defect formation is a slow process compared to ionic motion. On the other hand, a desorption mechanism becomes valid for negative applied voltages and reduces the velocity during the reset process. Since these processes act on different time scales but have to be taken into account within the concentrated model, a state-dependent activation energy is preferable. To this end, we introduce a state-dependent, normalized energy barrier

[TABLE]

with $\Phi_{\mathrm{a}}=\varphi_{\mathrm{a}}\,k_{\mathrm{B}}\,\vartheta$. The indices $\{0,1\}$ denote the boundary values corresponding to the state, whereas $\varphi_{\mathrm{a_{r}}}$ is the assumed activation energy for the reset process. The linear dependency with respect to the state variable can be adapted by other functions considering different material properties.

Beside the ordinary differential equation for the memristive behavior (2b), we have to formulate expressions for the algebraic input-output relation (2a) between voltage and current. The derivation of these relations is shown in the next subsection.

Input Output Relations

The DBMD consists of three state-dependent contiguous regions, see figure 8. The regions can be interpreted as single memristive systems, all of them coupled by the state variable. In order to classify this device as a memristive system, we have to describe the input-output relations for each region like in equation (2a). In the sequel, the individual regions are distinguished.

Schottky-contact

The Schottky-current can be described by the Schottky-equation. Hansen et al[6] has used an additional fitting parameter for the reverse Schottky-current to fit the measured data.

Taking a barrier lowering due to image charges [18] into account in combination with a state-dependent normalized Schottky-barrier height $\varphi_{\mathrm{s}}(z)$ and an ideality factor $n(z)$, the Schottky-current is given as

[TABLE]

Here, a dimensionless fitting parameter $\alpha_{\mathrm{f}}$ weights the barrier lowering term, which is caused by the Schottky-effect. The normalized Schottky-barrier thickness is

[TABLE]

where $d_{\mathrm{s}}$ denotes the absolute Schottky-barrier thickness and $D_{\mathrm{s}}$ is the normalization factor, with permittivity of vacuum $\epsilon_{0}$ and relative permittivity $\epsilon_{\mathrm{r}}$ of the electrolyte. The amplitude $I_{\mathrm{s}}=R_{\mathrm{i}}\,A\,\vartheta^{2}$ contains besides the temperature the effective Richardson constant $R_{\mathrm{i}}$ and the cross-sectional area $A$ of the device. Linear functions for the state-dependencies of the normalized Schottky-barrier height and the ideality factor have been assumed and can be expressed by

[TABLE]

with $\varphi_{\mathrm{s_{0}}}\leq\varphi_{\mathrm{s}}(z)\leq\varphi_{\mathrm{s_{0}}}$ and $n_{0}\leq n(z)\leq n_{1}$. The resulting input-output relation between current through and voltage drop over the Schottky barrier is given by equation (2k); it can be rewritten into the form of equation (2a) as

[TABLE]

thus defining resistance $\hat{R}_{\mathrm{s}}$. With this, we can identify the normal form of a memristive system in term of equations (2a) and (2b) for the Schottky-region.

Electrolyte region

Memristive behavior of the electrolyte region stems from the actual ion distribution within the electrolyte, which influences the conductivity [19]. A deeper modeling approach requires a state-dependent resistance for this region. For this purpose, we assume a linear state-dependency for the electrolyte region

[TABLE]

with $R_{\mathrm{e_{0}}}<\hat{R_{\mathrm{e}}}(z)<R_{\mathrm{e_{1}}}$ as the high and low resistance states, respectively. This in turn leads to the memristive normal form of equation (2a) and (2b) regarding the input-output relation

[TABLE]

Tunnel barrier

The last region to address is the Al2O3 layer. This layer is an electrically high quality tunnel barrier [6]. An electrical formulation of the tunnel current was given by Simmons [20]. There, a distinction depending on the applied voltage was introduced: low voltage regime $|u_{\mathrm{t}}|\approx 0$, intermediate regime $|u_{\mathrm{t}}|\leq\Phi_{\mathrm{t}}/q_{\mathrm{e}}$ and high voltage regime $|u_{\mathrm{t}}|>\Phi_{\mathrm{t}}/q_{\mathrm{e}}$, where $\Phi_{\mathrm{t}}$ denotes the tunnel barrier height and $u_{\mathrm{t}}$ the voltage drop over the tunnel barrier. Choosing the right equation for an accurate functionality of the concentrated model is important.

Kinetic Monte-Carlo simulations in [8] have shown that the restriction $u_{\mathrm{t}}<3~{}\mathrm{V}$ is fulfilled for the considered operating range $-2~{}\mathrm{V}\leq u\leq 3~{}\mathrm{V}$ of the device. This justifies the use of the intermediate Simmons equation for the concentrated model, because the tunnel barrier height of the DBMD is about $\Phi_{\mathrm{t}}\approx 3~{}\mathrm{eV}$, cf. [8]. A concentrated formulation of the Simmons equation leads to the tunnel current

[TABLE]

The normalized tunnel barrier height is denoted by $\varphi_{t_{0}}$. Here, $\alpha_{\mathrm{t}}$ is the normalized tunnel barrier thickness, where the absolute width $d_{\mathrm{t}}$ is normalized by $D_{\mathrm{t}}=\frac{\mathrm{h}}{4\pi\sqrt{2\,m_{\mathrm{e}}\,k_{\mathrm{B}}\,\vartheta}}$. The normalization constant is based on physical constants like Planck’s constant $h$ and mass of an electron $m_{\mathrm{e}}$.

Due to the Simmons formula [20], which accounts only for elastic tunneling, the effective barrier thickness depends on the ion distribution and with this on the internal state. This results in a state-dependent normalized tunnel barrier thickness. For the sake of an effective implementation, we have used a linear dependency, which leads to

[TABLE]

As it can be seen, the input-output relation of the tunnel barrier region reads

[TABLE]

which is a memristive system, cf. equations (2a) and (2b).

5 Simulation Results of the Concentrated Model

Several simulation models of memristive systems are available in the literature [17]. Most of them are based on mathematical descriptions, where the physical meanings of parameters are not obvious.

In this section, we want to verify the concentrated model by comparisons with kinetic Monte-Carlo simulations as well as measurements. To this end, an LTSpice implementation of the concentrated model is used. For investigations of small time scales $T\approx 100\,\mathrm{s}$ hysteresis curves are utilized. In contrast to that, long time-scales $T\approx 600\,\mathrm{s}$ are investigated by the step response of the device for two amplitudes. A semilogarithmic scale with respect to the absolute value of the current is used for the hysteresis curves. This improves the visibility of a gradual resistance change. To be consistent with measurements and kinetic Monte-Carlo simulations, the high resistance state $z=1$ is chosen as the initial state for all following simulations.

Short Time Scale Simulations: Hysteresis Curves

Short time scale investigations are done by considering resulting hysteresis curves for input voltages of the form depicted in figure 3 (inset). For the sake of consistency, the amplitudes are chosen as in measurements and kinetic Monte-Carlo simulations. The results are shown in figure 9, where for improved clarity the areas of some hysteresis loops have been shaded. In figure 9 a), hysteresis curves for exactly the same input voltage of figure 3 are shown. The good coincidence between LTSpice simulation, kinetic Monte-Carlo simulation and measurement is remarkable, regarding the complexity of the distributed model compared to the concentrated.

In figure 9 b)-d) hysteresis curves for different voltage amplitudes are depicted.

In coincidence with measurements and kinetic Monte-Carlo simulations, resulting curves for the concentrated model show that a threshold voltage $U_{\mathrm{T_{S}}}$ has to be reached to induce a memristive behavior, because for voltages below $U_{\mathrm{T_{S}}}$ no hysteresis occurs. To emphasize this behavior, figure 10 shows the voltage drops over individual regions with respect to the total voltage drop over the device.

As expected, until $U_{\mathrm{T_{S}}}$ is reached almost the whole amount of the applied voltage drops over the Schottky-region. When the voltage reaches $U_{T_{S}}$, the Schottky-contact becomes more and more a short-circuit. The forward direction of the diode leads to an increasing voltage drop over the electrolyte as well as over the tunnel barrier. The voltage drop over the electrolyte results in a change of the state variable and therefore in a decreasing of the device resistance, cf. [8].

Long Time Scale Simulations: Step Responses

For long time scale investigations, step responses of the device, with two different voltage amplitudes, have been measured and simulated. In figure 11, the step response for an input voltage amplitude of $2.5~{}\mathrm{V}$ (left) and $2.9~{}\mathrm{V}$ (right) is presented.

Again a good coincidence between concentrated and distributed model can be seen. This is comprehensible since the concentrated model is based on the distributed model. However, it should be stressed that further investigations, like the step response, were initiated by the concentrated approach. On top of that, the agreement with measured data is also important.

We see from simulation results, that the current is higher for an input voltage amplitude of $2.5~{}\mathrm{V}$ and lower for an input voltage amplitude of $2.9~{}\mathrm{V}$ compared to the measurement. Similar discrepancies can be observed in the hysteresis curve of figure 9 a). We have to mention that the model has of course an imprecision which could lead to such discrepancies.

Conclusion

In this work, we derived a concentrated model of a double barrier memristive device based on investigations of a distributed model and measurements. The concentrated model increases the mathematical accessibility, while maintaining physically meaningful parameters.

The DBMD has been classified as a memristive system, with the average ion position as its state, i. e. the drift velocity represents the memristive behavior. In this memristive system, a modified window function overcomes the boundary lock problem. For a gradual resistance change, physically justified state-dependencies of particular parameters have been utilized. Novel insights from actual research results, like adsorption and desorption mechanisms or ion formations within the electrolyte for long time scales, have been incorporated by an appropriate activation energy depending on both the state variable and the applied voltage. It is notable, that the presented approach can also be used for physical investigations. The step response of the DBMD was initiated during the modeling procedure of the concentrated model.

Although the distributed model is more complex than the concentrated one, LTSpice simulation results have shown a good coincidence compared with kinetic Monte-Carlo simulations as well as measurements. As an example, the simulation time with LTSpice was on the time scale of seconds, whereas a kinetic Monte-Carlo simulation with a distributed model takes approximately some hours.

This general approach in combination with physically meaningful parameters restrict the presented model not to a particular device. Instead, other resistive switching devices with different material compositions can be modeled only by adapting corresponding parameters. In conclusion, the concentrated model offers new possibilities for the investigation of complex neuromorphic circuits including real memristive devices, e.g. sensitivity analyses with respect to noisy conditions and parameter spread.

Acknowledgements

The financial support by the German Research Foundation (Deutsche Forschungsgemeinschaft - DFG) through FOR 2093 is gratefully acknowledged.