ABSTRACT (original article): The formation of metallic nanofilaments bridging two electrodes across an insulator is a mechanism for resistive switching. Examples of such phenomena include atomic synapses, which constitute a distinct class of memristive devices whose behavior is closely tied to the properties of the filament. Until recently, experimental investigation of the low-temperature regime and quantum transport effects has been limited. However, with growing interest in understanding the true impacts of the filament on device conductance, comprehending quantum effects has become crucial for quantum neuromorphic hardware. We discuss quantum transport resulting from filamentary switching in a narrow region where the continuous approximation of the contact is not valid, and only a few atoms are involved. In this scenario, the filament can be represented by a graph depicting the adjacency of atoms and the overlap between atomic orbitals. Using the theory of quantum graphs with locally diffusive node scattering, we calculate the scattering amplitude of charge carriers on this graph and explore the interplay between filamentary formation and quantum transport effects. CITATION (original article): Alison A. Silva, Fabiano M. Andrade, Francesco Caravelli (2024), Quantum graph models for transport in filamentary switching, arXiv:2404.06628. https://doi.org/10.48550/arXiv.2404.06628

Quantum graph models for transport in filamentary switching (QFCLaplace3)

This is the main function of our routine which connects all the other functions here defined.
The input parameters to evaluate the system are:
upot: the constant potential along the edges of the quantum graph;
Dv: the dimensions of the 3D space (Dv x Dv x Dv);
singleout: a parameter that allow us to switch between:
1: single output, only one system will be evaluated and save all of the data of this single iteration such as particles positions, electric field, and so on;
0: optimized to do many iterations and save only the main data as electrical current, voltage, etc., while discards the ones that will not be used.
Here we first define our electric potential difference function, which default is defined as
V
0
=0.35sin(0.5t)
, other parameter are defined such as the time step
(dt=0.1)
, number of points in the transmission probability evaluation
(Nk=600)
,...Then the initial system is defined with the function SetLeadsFlat3D, which defines both leads. With that, further calculations are made with the function PDens3DInter, where is evaluated the electric field and the filament composition by an assembly of atoms with overlapping bands. The next step is to transform this filament in a graph using GraphFromMat3, where the adjacent atoms in the system are used to define an adjacency matrix. From the adjacency matrix, we can obtain the scattering transmission in a quantum graph defined by it, for this purpose the function TQG is used. From that is possible to evaluate the current from the Landauer-Büttiker formalism, from which we got the formula for the system with two scattering channels
I(V)
G
0
∞
∫
-∞
dE(
f
L
(E)/N-
f
R
(E)/M)
h
2
2
e
σ
0
(V)+T(E)
where
G
0
is the conductance quantum,
f
L
(E)
and
f
R
(E)
are obtained from the Fermi-Dirac distribution, which can be evaluated with the chempot function,
σ
0
is an additional conductivity of the device, which can be either constant or voltage-dependent, and
T(E)
is the transmission probability obtained from the quantum graph scattering transmission evaluation. From these evaluations, the function returns Vt: a list containing the voltage for each step of time, It and Ifp: the current as function of time for a trivial or a Frenkel-Poole insulator, respectively, and Ampt: a list with the scattering transmission probability of the quantum graph modeled by the filament for each step of time.

Initial boundaries (SetLeadsFlat3D)

In this function the objective is to define the initial boundaries of the system without filaments. Here we first define the dimensions of the system with the parameters H, W and Dv, and then we define the top and bottom planes of the cube as two plates which are obtained by the following function:
In[]:=
SetLeadsFlat3D[H_,W_,Dv_]:=Module[{NpTop=ConstantArray[0,{W,Dv,H}],NpBottom=ConstantArray[0,{W,Dv,H}]},​​NpTop〚All,All,1〛=ConstantArray[1,{H,W}];​​NpBottom〚All,All,-1〛=ConstantArray[1,{H,W}];​​{NpTop,NpBottom}];
It will return a list of two 3D arrays, one for each plane, which will define the leads of our system.

Example:

Lets define a system with dimensions 20 x 20 x 20, thus we have
In[]:=
exSetLeadsFlat3D=SetLeadsFlat3D[20,20,20];
Now can we visualize them with
In[]:=
GraphicsRow[ArrayPlot3D/@exSetLeadsFlat3D]
Out[]=

Electric Field evaluation (LaplaceSolver3D)

In this function we aim to obtain the electric field between the two planes that define the leads in our system. To obtain that we need the following parameters:
V0: the electric difference potential between the leads;
MyBoundaryTop and MyBoundaryBottom: the two leads with +V0/2 and -V0/2 electric potentials respectively;
First we define a 3D array containing the initial potential in all of the discrete points of the system, being
+
V
0
/2
in the top plane,
-
V
0
/2
in the bottom one and
0
elsewhere. Then for all points, except the ones in the border of the defined space, we evaluate its electric potential by taking the mean value of six neighbor points, two for each direction (one before and one after). This iteration will be realized until the difference between the new and the last iteration is less or equal than the defined error parameter.Finally the electric field is evaluated from the potential matrix, by taking the differences between the electric potential values in each directions, as an discrete gradient. The function constructed for this purpose is the one as follows:
In[]:=
LaplaceSolver3D[V0_,MyBoundaryTop_,MyBoundaryBottom_]:=Module{xdim,ydim,zdim,VNow,VPrev,error,ex,ey,ez},​​xdim=Dimensions[MyBoundaryTop]1-1;​​ydim=Dimensions[MyBoundaryTop]2-1;​​zdim=Dimensions[MyBoundaryTop]3-1;​​VPrev=ConstantArray[0,{xdim+1,ydim+1,zdim+1}];​​VNow=
1
2
V0(MyBoundaryTop-MyBoundaryBottom);​​error=Max[Abs[VNow-VPrev]];​​VPrev=VNow;​​Whileerror>0.01`,​​Do​​IfMyBoundaryTopi,j,k0&&MyBoundaryBottomi,j,k0,VNowi,j,k=
1
6
(VPrev〚i-1,j,k〛+VPrev〚i+1,j,k〛+VPrev〚i,j-1,k〛+VPrev〚i,j+1,k〛+VPrev〚i,j,k-1〛+VPrev〚i,j,k+1〛);;​​,{i,2,xdim},{j,2,ydim},{k,2,zdim};​​error=Max[Abs[VNow-VPrev]];​​VPrev=VNow;;​​MatrixPosaux=VNow;​​ez=MapAt[0.5`#1&,ArrayPad[Differences[VNow,{0,0,1}],{{0,0},{0,0},{0,1}}]+ArrayPad[Differences[VNow,{0,0,1}],{{0,0},{0,0},{1,0}}],{All,All,2;;-2}];​​ey=MapAt[0.5`#1&,ArrayPad[Differences[VNow,{0,1,0}],{{0,0},{0,1},{0,0}}]+ArrayPad[Differences[VNow,{0,1,0}],{{0,0},{1,0},{0,0}}],{All,2;;-2,All}];​​ex=MapAt[0.5`#1&,ArrayPad[Differences[VNow,{1,0,0}],{{0,1},{0,0},{0,0}}]+ArrayPad[Differences[VNow,{1,0,0}],{{1,0},{0,0},{0,0}}],{2;;-2,All,All}];​​{ex,ey,ez};
The electric field is returned in the form {ex, ey, ez} being a list of three 3D arrays, each of these arrays the component for the electric field in the directions x, y and z for each discrete point of the defined system.

Example:

Lets define a
20
x
20
x
20
system by using the function SetLeadsFlat3D, which will give us the two leads.
In[]:=
exSetLeadsFlat3D=SetLeadsFlat3D[20,20,20];
Also we will plot both leads as top and bottom planes, with potential
+0.5
and
-0.5
, respectively.
In[]:=
ArrayPlot3D[Transpose[.5exSetLeadsFlat3D[[1]]-.5exSetLeadsFlat3D[[2]],1<->3],PlotLegends->Automatic]
Out[]=
-0.5
0.5
Now let us evaluate the electric field between these leads with the LaplaceSolver3D, by defining V0=1:
It is interesting to rearrange the outputs in a single array,
this allow us to plot the electric field with the following command :

Particles interaction and filaments construction (PDens3DInter)

This function aims to model the filaments construction, for that we must define the initial parameters:
treetopl and treebottoml: list of coordinates {x,y,z} for each component of the filament, which will become a vertex of the graph; Initially both are defined as an empty list and will grow as this function is re utilized in a loop; The difference between top and bottom are due the particle will start from different leads according to the electric potential in each lead;
MatrixPos: An array with the dimension of all the possible positions. As a particle become part of the filament, an label will be give to it and stored in this array in the same position as the particle;
nnodes: number of vertices in the graph that model the filament;
T: a time step parameter;
H and W: parameters related with the dimension of the system, W x W x H (along with this routine all of three are the same);
pdecay: parameter related to the probability of an vertex be removed of the graph or the decay of some part of the filament;
V: the electric potential difference at that step of time;
Ppos: Position of all the particles that are not connected to the filament (free particles between the leads), this list contains for each particle, its early position and its actual position;
TreeTop0, TreeBottom0, TreeTop and TreeBottom: the early and new matrices which storage which discrete positions are already occupied and are used to verify if the free particles are close enough to be attached to the filament.

Graph construction from filament position(GraphFromMat3)

Transmission in a Quantum Graph (TQG)

Through this method, the Green function of this system is obtained as

Example:

First let us define a graph like the Wheel graph with 5 vertices:
Now we take its adjacency matrix and print it:

Other quantities and functions (qe, me, ℏ, chempot, DecayM2)

Examples and data plots

Single iteration with video generation

First we define the list with the data that we will use and run the function QFCLaplace3 with the target parameters, e.g. here will be defined as potential in the edges upot=0, Dv=20 and singleout=1 to save all the data for a single iteration, which includes electric potential vectors, particles positions,... All the data is automatic stored in the same file as this notebook.
After running this, now we can generate the videos from these data. Following we will describe each of the plots and the code to generate them.
Note: “MaTeX” package is needed.
Voltage versus time
Current versus time
Current (Frenkel-Poole insulator) versus time
Current versus voltage
Current (Frenkel-Poole) versus time
Transmission
Filament formation
Plots combine

CITE THIS NOTEBOOK

Quantum graph models for transport in filamentary switching​
by Alison A. Silva, Fabiano M. Andrade & Francesco Caravelli
Wolfram Community, STAFF PICKS, December 3, 2024
​https://community.wolfram.com/groups/-/m/t/3332704