# Composition-driven Mott transition within SrTi1-xVxO3

A. D. N. James  H. H. Wills Physics Laboratory, University of Bristol, Tyndall
Avenue, Bristol, BS8 1TL, United Kingdom  M. Aichhorn  Institute of
Theoretical and Computational Physics, TU Graz, NAWI Graz, Petersgasse 16,
8010 Graz, Austria  J. Laverock  School of Chemistry, University of Bristol,
Cantocks Close, Bristol, BS8 1TS, United Kingdom

###### Abstract

The last few decades has seen the rapid growth of interest in the bulk
perovskite-type transition metal oxides SrVO3 and SrTiO3. The electronic
configuration of these perovskites differs by one electron associated to the
transition metal species which gives rise to the drastically different
electronic properties. Therefore, it is natural to look into how the
electronic structure transitions between these bulk structures by using
doping. Measurements of the substitutional doped SrTi1-xVxO3 shows an metal-
insulator transition (MIT) as a function of doping. By using supercell density
functional theory with dynamical mean field theory (DFT+DMFT), we show that
the MIT is indeed the result of the combination of local electron correlation
effects (Mott physics) within the t2g orbitals and the atomic site
configuration of the transition metals which may indicate dependence on site
disorder. SrTi1-xVxO3 may be an ideal candidate for benchmarking cutting-edge
Mott-Anderson models of real systems. We show that applying an effective
external perturbation on SrTi1-xVxO3 can switch the system between the
insulating and metallic phase, meaning this is a bulk system with the
potential use in Mott electronic devices.

##  I Introduction

There is continuing interest in oxide based materials as they exhibit exotic
phenomena, such as colossal magnetoresistance; high-temperature
superconductivity; and Mott insulating phases, or respond to external stimuli
which is highly desirable for innovating electronic devices [1, 2, 3]. In
particular, perovskite-type transition metal oxides ABO3 show many desirable
properties such as a metal-insulator transition (MIT) [4] which are of
interest for developing Mott devices [5, 3]. It is the many-body correlated
behaviour of electrons in these materials which gives rise to the emergence of
these exotic phenomena. Here, the interaction energy between electrons is
comparable to their kinetic energy. Therefore, it is important to understand
the correlated behaviour of electrons and how to manipulate them for future
technological applications.

Bulk cubic SrVO3 is a prototypical correlated material [6] and is thought to
have several applications owing to its high performance [7, 8, 3]. It has a
three peak spectral structure around the Fermi level with the central part
consisting of sharp $3d$ quasiparticle bands at low excitation energies giving
a well-defined Fermi surface [9, 10]. These quasiparticle bands are
accompanied by localised states (the so called incoherent Hubbard sidebands)
at an energy scale comparable with the Coulomb repulsion parameter, $U$ [11,
12, 13]. This material has been the subject of many density functional theory
with dynamical mean field theory (DFT+DMFT) studies [6, 11, 14, 15, 16, 17],
where only one electron exists within the V $3d$ $t_{\rm 2g}$. Often only the
V $3d$ $t_{\rm 2g}$ states around the Fermi level are used in the DMFT
calculations as the $e_{\rm g}$ do not hybridise with the $t_{\rm 2g}$ states
and are also unoccupied. DFT+DMFT is able to describe all of the on-site local
correlations [18, 15], and has been well-tested on SrVO3 with very good
results, including the energetics and spectral weight of Hubbard sidebands and
QP renormalisation [11, 14, 15, 16, 17]. However, $GW$+DMFT [19] predictions
reinterprets SrVO3 as a weakly correlated material with low static local
interactions, since their results show pronounced plasmonic satellites due to
screening.

The importance of electron correlations in a system may be gauged in terms of
the ratio, $U/W$, where $W$ is the non-interacting bandwidth as predicted by
DFT. There are certain factors which been shown to have influenced this ratio.
Crystal field effects have been a prime example of this as it can reduce $W$
to give the Mott insulating phase as predicted by DFT+DMFT of bilayer capped
SrVO3 [20] which agrees with that seen in experimental spectral data of SrVO3
thin films (below approximately 6 unit cells) [21, 22]. These effects are also
predicted in other vanadates where these can be tuned by strain [23, 24, 25,
26, 27]. Conversely, $W$ can also be influenced by quantum confinement of
heterostructure superlattices [28] in agreement with the corresponding
experimental data [29, 30, 31, 32].

Bulk SrTiO3 is a band-insulator despite sharing the same cubic structure as
SrVO3 at room temperature (below $\sim$ 105 K, SrTiO3 exhibits a tetragonal
structural transition [33]). SrTiO3 has been the focus of modern day research
due to its many desirable properties which can be exploited for developing
electronic devices [1, 34, 35, 36, 37]. The main interest in SrTiO3 is it
often being used as a key building block for new oxide heterostructures either
as a substrate to grow other oxide compounds on or within oxide
heterostuctures which leads to unique phenomena originating from the interface
between oxides (examples of which are discussed in Refs. [38, 39, 40]).
Evidently, it is this versatility of SrTiO3 which has resulted in its own
field of research.

It is interesting how these two bulk (structurally simple) perovskite
transition metal oxides compounds where the transition metal species
effectively differs by an electron (and a proton) give rise to such different
electronic structures. Therefore, tuning the electronic structure between
SrVO3 and SrTiO3, achieved by introducing impurities or defects in the
material, will give an insight into the rise of interesting correlated
electron behaviour. For example, substitutional doping has been shown to be an
effective way to introduce filling-controlled MIT in the aliovalent A-site
substitution of the intensely studied La1-xSrxVO3 [41, 42] (where $\rm x$
represents doping concentration). Recent studies of B-site substitutionally
doped SrTi1-xVxO3 show that there is a temperature-dependent MIT which is
composition driven [43, 44, 45, 46]. Beyond this composition-driven MIT, other
properties of this doped system have been investigated such as its impedance
spectroscopy and electrical conduction mechanism for certain dopings [47] and
it may also be classified as a correlated transparent conductive oxide [44]
which has the potential for use in optoelectronic devices. Several studies
have measured the critical doping for the MIT to be at different dopings,
ranging from between 0.4 [44], 0.6 [46], 0.67 [45] and 0.7 [43]. The critical
doping was shown to be restrained between 0.5 and 1 when introducing oxygen
vacancies [48].

The driving force behind the MIT in SrTi1-xVxO3 is still not fully understood
with speculation of its mechanism. Tsuiki et al. interpreted their findings by
arguing that the distortion from the Jahn-Teller effect [49] makes the donor
levels so deep that the vanadium $d$-orbitals cannot overlap to form a band
[46]. Hong et al. proposed that this is MIT is controlled by the Anderson
localized states in which a MIT occurs where the mobility edge crosses the
Fermi level [43]. Subsequent studies then proposed that the mechanism behind
this composition drived MIT is due to a combination of both electron
correlation effects (Mott physics) and disorder-induced (Anderson)
localisation [44, 45], where the Ti4+ ion substitution introduces Anderson-
localized states along with lattice distortions which result in a reduction in
the effective bandwidth $W$. Therefore, it is clearly necessary to understand
how the electronic structure transitions in SrTi1-xVxO3 between the band-
insulator to the correlated metal from a theoretical perspective in order for
this doped system to be effectively manipulated.

Here, we look at the origin of the dominant mechanism which gives rise to the
MIT in SrTi1-xVxO3. We focus on using the DFT+DMFT method on supercells of
different configurations of the highly symmetric atomic positions in order to
gain insight into the mechanism. We find that from the set of configurations
studied here, the MIT is highly dependent on both the local environment of the
V atoms in the doped systems (which is influenced by the site configuration)
and the Mott local electron correlation effects within the t2g orbitals. The
doping-dependent non-interacting $W$ in combination with the on-site $U$
pushes the system, depending on the $U$/$W$ ratio, into either the metallic
(low $U$/$W$) or insulating phase (high $U$/$W$). Around the MIT critical
doping, the MIT appears to be sensitive to crystal field effects on $W$
introduced from site configuration of the transition metal species within the
unit cells. Finally, we show that the effect of an external perturbation on
Mott insulating dopings can push them metallic which has potential
applications for Mott devices.

##  II Method

Figure 1: The unit cells of SrTi1-xVxO3 used in the supercell DFT and DFT+DMFT
calculations. The configurations and stoichometries of SrTi1-xVxO3 are shown
in (a) – (d). The labels for each species are shown in (b). Bulk SrTiO3 has
the same crystal structure as SrVO3 shown in (a).

DFT calculations were performed with the full potential augmented plane-wave
plus local orbitals (APW+lo) elk code [50], using SrTi1-xVxO3 cubic lattices
shown in Fig. 1 as input. These configurations were chosen to ensure inversion
symmetry and maximise V-V and Ti-Ti neighbour distances within the $2\times
2\times 2$ cell (other structures were also considered and are discussed in
Appendix A and during the text). The DFT calculation used the Perdew-Burke-
Ernzerhof (PBE) generalized gradient approximation (GGA) for the exchange-
correlation functional [51] and was converged on a $12\times 12\times 12$
Monkhorst-Pack $\mathbf{k}$-mesh of 84 irreducible $\mathbf{k}$-points in the
irreducible Brillouin zone. The (simple cubic) lattice parameter $a$ at which
the DFT total energy is minimised for each of these structures was used. For
comparison, additional virtual crystal approximation (VCA) calculations were
performed at the same doping levels, ${\rm x}$, by replacing the V ion in a
single SrVO3 unit cell with an effective Ti1-xVx ion. We emphasise that VCA
works well to describe compositional alloying and disorder of atoms which do
not contribute to the states at the Fermi energy, and is not expected to
perform well here.

Figure 2: The DFT t2g partial density of states (PDOS) within the transition
metal muffin-tin for each of the supercells shown in Fig. 1, including bulk
(a) SrVO3 and (e) SrTiO3. At the DFT level, the Fermi level is located within
the PDOS of the V t2g orbitals for the systems contains V, whereas it is
within the band-gap for SrTiO3.

For the DFT+DMFT calculations, the elk code [50] was used in combination with
the toolbox for research on interacting quantum systems (TRIQS) library [52].
This so-called elk+TRIQS DFT+DMFT framework is described in Ref. [53]. The DFT
outputs were interfaced to the TRIQS/DFTTools application of the TRIQS library
[54] by constructing Wannier projectors, as described in Ref. [53], derived
from the bands around the Fermi level which are predominantly V t2g in
character (see Fig. 2). The V t2g Wannier charge of each structure is 1
electron per V site. The DMFT part of the DFT+DMFT calculation was implemented
using the continuous-time quantum Monte Carlo (CT-QMC) solver within the
TRIQS/CTHYB application [55] with the Kanamori-Hubbard interaction Hamiltonian
parameterised by the Hubbard interaction $U=4.0$ eV and Hund exchange
interaction $J=0.65$ eV, unless otherwise specified. We approximated the
double counting within the fully localised limit. These $U$, $J$ values and
double counting approximation are similar to those used for previous bulk
SrVO3 DFT+DMFT calculations [56]. We used the fully-charge-self-consistent
(FCSC) DFT+DMFT method with a total of $1.12\times 10^{8}$ Monte Carlo sweeps
within the impurity solver for each DMFT cycle. Owing to symmetry, there is
only one impurity to solve for each doped system shown in Fig. 1. An inverse
temperature $\beta=40$ eV-1 ($\sim$ 290 K) was used. The spectral functions
were calculated by analytically continuing the DMFT self-energy obtained from
the LineFitAnalyzer technique of the maximum entropy analytic continuation
method implemented within the TRIQS/Maxent application [57].

In order to characterise the metal-insulator transition, certain quantities
were extracted from the self-energies $\Sigma$ and Green’s functions, which we
define here for clarity. The presented spectral functions at the Fermi level
are determined directly from the imaginary time Green’s function by

$A(\omega=0)=\frac{\beta G(\tau=\frac{1}{2}\beta)}{\pi},$ (1)

where $\beta$ is the inverse temperature in natural units. The $A(\omega=0)$
here is an averaged quantity over a frequency window approximately equal to
$\beta^{-1}$. These are normalised to the bulk SrVO3 spectral function at the
Fermi level to give $\overline{\rm A}$($\omega=0$). Our calculations (here and
in Appendix C) where we refer to applying an effective potential (gate voltage
Vg) were achieved by providing a fractional excess electron charge to the
$\rm{x}=0.5$ (and 0.25) supercell, resulting in a chemical potential shift of
the uncorrelated states of eVg.

##  III Results

Figure 3: The DFT+DMFT integrated V Wannier spectral function A($\omega$)
(including the V-Ti hybridised states in the Wannier projectors) of the
supercells with decreasing V concentration from (a) to (d). The metal-
insulator transition (MIT) occurs close to ${\rm x}\sim 0.5$. The Fermi level
has been placed at the centre of the band-gap for the insulating systems.

Figure 4: The DFT+DMFT k-resolved spectral function A(k, $\omega$) of
SrTi1-xVxO3, with decreasing V concentration from (a) to (d). Solid lines
indicate the non-interacting DFT band structure. The bulk A(k, $\omega$) and
band structure in (a) have been back-folded into the $2\times 2\times 2$
supercell Brillouin zone.

The non-interacting t2g partial densities of states (PDOS) from DFT of
SrTi1-xVxO3 are shown in Fig. 2. Except for SrTiO3, which has no occupied
electrons in the $t_{\rm 2g}$ manifold, all V-doped systems are metallic
within DFT. The V t2g PDOS of the supercells in Figs. 2 (b)–(d) clearly lie in
an energy range near $E_{\rm F}$, which we refer to as the non-interacting
bandwidth. The non-interacting V $t_{\rm 2g}$ bandwidth decreases as the V
concentration $\rm x$ is decreased, dropping by a factor of 3 between x = 1
and x = 0.25. There are additional contributions in the Ti and V t2g PDOS
which arise due to hybridisation between these orbitals. We note that in the
DFT calculations the V t2g PDOS is the largest contribution to the states at
and around the Fermi level in all the systems and configurations investigated.
Although DFT adequately qualitatively describes the quasiparticle peak of the
V t2g states in SrVO3, it fails to capture the quasiparticle bandwidth and the
Hubbard bands due to its neglect of electron-electron interactions, and no MIT
is observed at the level of DFT for the doped systems. We note that recent
DFT+U calculations of SrTi0.75V0.25O3 in the low-temperature (orthorhombic)
structure of SrTiO3 were shown to be insulating due to distortion-induced
crystal-field effects [58].

Figure 3 shows the results of supercell DFT+DMFT calculations, where the
Wannier V DOS is shown for each doped system. In agreement with previous
calculations [6, 56] the DFT+DMFT correctly predicts the Hubbard bands for
pure SrVO3 (Fig. 4), which appear at $\approx-1.7$ eV (lower Hubbard band,
LHB) and $\approx 2.7$ eV (upper Hubbard band, UHB). As the Ti concentration
is increased, the V Hubbard band intensity rises as spectral weight is
transferred from the quasiparticle peak, and the bandwidths of the V
quasiparticle states are renormalised with respect to the non-interacting
(DFT) bandwidth (see Fig. 4). At a critical doping, $\rm{x}_{\rm c}\sim 0.5$,
the system undergoes a metal-insulator transition as the spectral weight is
completely diminished at the Fermi level. This is in very good agreement with
previous experiments. Similar results are obtained with a different value of
$U$ and for different atomic arrangements of the supercell (see Appendix A and
following discussion). We note that the atomic configuration of the
$\rm{x}=0.75$ supercell employed here lifts the degeneracy of the V t2g
orbitals, which gives rise to a small dip in the V DOS in Fig. 3(b); this is
an artefact of the supercell approach. Amongst the upper Hubbard bands of the
doped systems, there are additional distinct broad features above about 2.5
eV. These correspond to hybridised V t2g and Ti t2g states as seen in the DFT
PDOS in Fig. 2.

Figure 5: (a) The VCA and supercell DFT (normalised) orbital-averaged
bandwidths $\overline{W}$ of the V t2g bands located around the Fermi level
(which relate to the bandwidths of the V t2g PDOS shown in Fig. 2). (b) The
VCA and supercell DFT+DMFT (normalised) integrated V Wannier spectral function
at the Fermi level $\overline{\rm A}$($\omega=0$) showing the metal-insulator
transition near $\rm{x}=0.5$. The overline over the variables shown in this
figure refers to the quantity being normalised to the corresponding bulk SrVO3
value. These quantities highlight the reduced coherence of the V t2g states
near the transition close to $\rm{x}=0.5$. The red (left) to blue (right)
background colour of (b) is a guide to the eye of which regions are insulating
(left) and metallic (right).

As illustrated in Figs. 4 and 5 (a), all doped structures are metallic from
the DFT predictions, which clearly shows that (near $\rm{x}_{\rm c}$) the
insulating phase is Mott in nature, arising from the local electron
correlations included in DMFT. The supercell non-interacting (orbital-
averaged) bandwidth from bare DFT ($\overline{W}$) is shown in Fig. 5 (a), and
decreases almost linearly from SrVO3 to SrTi0.75V0.25O3. This reduction in the
quasiparticle bandwidth is greater than that achieved in other doped SrVO3
systems, such as SrxCa1-xVO3 [59, 13]. Correspondingly, the integrated
spectral function at the Fermi level $\overline{\rm A}$($\omega=0$) [displayed
in Fig. 5 (b)] increases as the quasiparticle peak width shrinks as the doping
approaches $\rm{x}_{\rm c}$, indicating the increasing influence of the local
electron correlations within the t2g orbitals (and has been seen for tuning
$U$ through the MIT in SrVO3 superlattice heterostructures [28]). These can be
explained by the $U$/$W$ ratio increasing when approaching the MIT from high
doping in the metallic phase. In the insulating phase, $\overline{\rm
A}$($\omega=0$) is equal to zero as expected. Figs. 5 (a) and (b) show our VCA
DFT and DFT+DMFT calculations used as an alternative to approximate the effect
of the substitutional doping. However, Fig. 5 (b) shows that these
calculations do not give rise to an MIT. Here, the effect of the VCA is to
effectively shift the t2g states (with respect to the Fermi level) to conserve
charge, and hence the almost constant $\overline{W}$ with respect to doping in
Fig. 5 (a). The VCA unsurprisingly does not capture the variations in the
potentials around each transition metal site which evidently needs to be
included.

As with the results of monolayer [23], bilayer [20] or superlattice
heterostructures [28] of SrVO3, the spectral weight from the V t2g states in
the $\rm{x}=0.25$ and $0.5$ doped systems appears to fully redistribute to the
Hubbard bands in this insulating phase. The UHB appears to have a momentum
dependent dispersion in Fig. 4 (c) and (d), similar to that seen in Ref. [28].
This dispersive UHB, which is strongly reminiscent of the underlying QP
dispersion, is very prominent in the A(k, $\omega$) of $\rm{x}=0.25$ system.
DFT+DMFT calculations using a real-time impurity solver have shown that the
UHB of bulk SrVO3 is composed of three atomic states (namely the excited two
electron [$N=2$] atomic states) [60]. The ill-posed nature of analytical
continuation results in the finer details of these atomic states being “washed
out” to give the typical broad UHB feature seen in bulk SrVO3 [60] (see Fig. 3
(a)). In our results, these atomic states appear to be more prominent in the
spectral function derived from the analytically-continued self-energy for the
low $\rm{x}$ doped systems (also see Appendix B). Indeed, the shouldered
feature around 1 eV in the $\rm{x}=0.25$ system appears to be composed of the
two singly occupied $N=2$ atomic states (see Appendix B).

Figure 6: The (a) 1$\times$1$\times$2 and (b) 2$\times$2$\times$1 supercell
configurations investigated for $\rm{x}=0.5$ doping. The labels for each
species are shown in (b). (c) $\overline{\rm A}$($\omega=0$) (as defined in
Fig. 5) for different V and Ti site configurations within different unit cells
of SrTi0.5V0.5O3, given as a function of the corresponding normalised DFT
bandwidth $\overline{W}$. The red (left) to blue (right) background colour of
(c) is a guide to the eye of which regions are insulating (left) and metallic
(right).

The impact of atomic configuration is clearly shown for $\rm{x}=0.5$ in Fig.
6. In this figure, three different unit cells of $\rm{x}=0.5$ are presented
with the V atoms having the shortest average intra-species distance in the
1$\times$1$\times$2 unit cell (Fig. 6 (a)), the longest average intra-species
distance in the 2$\times$2$\times$2 (i.e., the structure in Fig. 1 (c)) and
the average intra-species distance in the 2$\times$2$\times$1 unit cell [Fig.
6 (b)] is inbetween the other two. The corresponding results of these unit
cells presented in Fig. 6 (c) show that the $\rm{x}=0.5$ is at the precipice
of the MIT and is highly dependent on the position configurations of the
transition metal species. It also clearly indicates the critical normalised
non-interacting orbital-averaged bandwidth value $\overline{W}_{\rm c}$ (which
is around 0.55 eV) at which the MIT occurs. All of the $\rm{x}=0.5$
configurations can be driven into the insulating phase with a modest increase
in $U$ to 5 eV (while $\rm{x}=0.75$ remains metallic). Experimentally, the
rather broad range of $\rm{x}_{\rm c}$ that is reported (between 0.4 and 0.7)
is also influenced by other factors beyond disorder such as the different
substrates (e.g. SrTiO3 [44] or LSAT [45]) used, film thicknesses, and
measuring bulk crystals [46, 43]. In particular, SrTiO3 and LSAT substrates
introduce quite different in-plane strain to the SrTi1-xVxO3 lattice in
epitaxial thin films, so that crystal-field effects become important near
$\rm{x}_{\rm c}$.

Figure 7: The Wannier V t2g A($\omega$) with the V-Ti hybridised states of the
$\rm{x}=0.5$ system [in Fig. 3 (c)] and of the same system but doped with an
additional 0.025 of an electron per V site (e-/site). This is equivalent to a
gating voltage Vg of 0.26 V. This doping pushes the system into being metallic
showing that external perturbations with the same effect can switch the
$\rm{x}=0.5$ system between the metallic and insulator state. This gives the
possibility of SrTi1-xVxO3 being used for Mott electronic devices. The Fermi
level has been placed at the centre of the band-gap for the insulating system.

Finally, we show that these substitutionally-doped SrTi1-xVxO3 system near
$\rm{x}_{\rm c}$ can be switched between the metallic and insulating state in
Fig. 7. By applying an external (gating) potential Vg (which gives the
equivalent effect as the charge excess in these charge excess calculations) in
the $\rm{x}=0.5$ system, this results in it transitioning from the Mott
insulating phase to the metallic one. The Vg for $\rm{x}=0.5$ is similar to
that in the bilayer system presented by Zhong et al. [20]. A larger potential
is needed for the $\rm{x}=0.25$ system (discussed in Appendix C) as expected
as it is further away from $\rm{x}_{\rm c}$. This shows that these systems
could be used in Mott devices with this switching mechanism controlled by an
appropriate external perturbation. As these substitutionally-doped systems are
bulk in nature, these will not be prone to issues associated with maintaining
a pristine surface as would be required for the bilayer-capped SrVO3 system
presented by Zhong et al. [20].

Given our DFT+DMFT supercell results, it is clear that more detailed
investigations near $\rm{x}_{\rm c}$ will be invaluable. Computationally,
calculations of doping levels at a finer granularity than presented here
remain out of reach using an all-electron supercell approach. However, new
methods beyond DFT+DMFT that are appropriate for disordered systems are
rapidly becoming available. For example, Refs. [61, 62] describe new
approaches that are able to capture the local atomic configuration effects as
well as incorporating the additional Anderson localisation effects due to
site-disorder for interacting electron systems. Features in the spectral
function such as localised Anderson states (originating from disorder
partitioning the electron states into either localised or extended states)
within the Mott insulating gap are present in such models (e.g., see Ref.
[61]) which are not captured in the supercell results presented here. With the
clear dependence on the site configurations and the local electron correlation
effects of the V t2g states around the Fermi level, SrTi1-xVxO3 could serve as
a prototypical test material for initial real-material calculations of solving
such Anderson-Mott systems, just as SrVO3 was an initial focus when developing
DFT+DMFT.

##  IV conclusion

We have shown that the measured trends around the MIT in the doped SrTi1-xVxO3
system are captured well by DFT+DMFT predictions of different supercell
configurations. Refinements in both DFT and DMFT implementations allows these
supercell DFT+DMFT calculations to be more computational feasible since when
some of these ARPES measurements were published about a decade ago. From our
supercell results, it is clear that both the local electron correlations and
site configurations are vital for this MIT. It is the site configurations
which influences the non-interacting bandwidth, which in combination with the
electron correlation effects from the onsite Coulomb interaction gives rise to
the MIT. This therefore is a likely candidate for a Mott-Anderson MIT. We also
note that atomic states associated to the UHB are more prominent in these
doped systems, whereas in bulk SrVO3 these atomic states are “washed out” as a
consequence of analytic continuation, yielding the familiar broad UHB feature.

By providing an external perturbation, equivalent to the effect of an
additional small excess of electron doping, SrTi1-xVxO3 can be made to switch
between the Mott insulating and metallic phases, giving another possible
avenue to develop Mott electronic devices. These results scratch the surface
of the complex strong local correlated electron behaviour meaning that further
work is needed to probe this possible Mott-Anderson MIT to understand how the
different factors such as site disorder and the local electron correlations
influences the MIT around the critical doping which is needed to realise its
potential for Mott devices.

##  V Acknowledgements

A.D.N.J. acknowledges the Doctoral Prize Fellowship funding and support from
the Engineering and Physical Sciences Research Council (EPSRC). Calculations
were performed using the computational facilities of the Advanced Computing
Research Centre, University of Bristol (http://bris.ac.uk/acrc/). The VESTA
package (https://jp-minerals.org/vesta/en/) has been used in the preparation
of some figures.

##  Appendix A Unit cell configuration

Figure A.1 shows all the unit cell configurations (up to 2$\times$2$\times$2)
used for this work. The lattice parameter was determined by minimising the
total DFT energy as discussed in the main text. For the $U=4.0$ eV DFT+DMFT
calculations, each configuration for the $\rm{x}=0.25$ doped systems are all
insulating whereas all of the $\rm{x}=0.75$ (and 1.0) configurations are
metallic. For the $\rm{x}=0.5$ configurations, we see that the
metallic/insulating phase is dependent on the configuration (which influences
the DFT “non-interacting” bandwidths $W$ and as such on the electron
correlations which are usually quantified by the $U$/$W$ ratio $U$/$W$ ratio)
as discussed in the main text. We find that using $U=5.0$ eV pushes all of the
$\rm{x}=0.5$ configurations insulating as a consequence of its effect on the
electron correlations seen by the $U$/$W$ ratio. We note that the band gap in
the insulating phase increases as the $W$ decreases, again as a result of
enhancement of local correlation quantified by the $U$/$W$ ratio.

Figure A.1: The different investigated SrTi1-xVxO3 configurations up to 2x2x2
sized unit cells for $\rm{x}=0.25$, 0.5 and 0.75. The SrVO3 unit cell is
included for completeness. Each configuration is masked with an translucent
box which indicates whether the DFT+DMFT $U=4.0$ eV results of that
configuration is metallic (light blue) or insulating (light red), for example,
the SrVO3 unit cell is metallic.

##  Appendix B $U$-varied Metal-Insulator Transition

Figure B.1: The effect of varying $U$ (for fixed $J=0.65$ eV) on the
$\rm{x}=0.25$ Wannier V t2g spectral functions (calculated from the Wannier
projectors without the Ti-V hybridised states) through the metal-insulator
transition. The Fermi level has been placed in the centre of the band-gap for
the insulating results.

We also show the metal-insulator transition (MIT) from varying the on-site
Hubbard $U$ (for fixed $J=0.65$ eV) in Figs. B.1 and B.2 for doping
$\rm{x}=0.25$ and $0.5$ (of the unit cell configuration in the main text),
respectively. We see from the Wannier spectral functions in Figs. B.1 and B.2
(which are calculated from using the Wannier projectors constructed from the V
t2g bands at the Fermi level) that the spectral weight of the quasiparticle
states is totally depleted and redistributed to the upper and lower Hubbard
bands for both the $\rm{x}=0.25$ and $0.5$ doped systems. The Mott insulating
band gap increases with $U$ as expected for these Hubbard bands.

Figure B.2: The effect of varying $U$ (for fixed $J=0.65$ eV) on the
$\rm{x}=0.5$ Wannier V t2g spectral functions (calculated from the Wannier
projectors without the Ti-V hybridised states) through the metal-insulator
transition. The Fermi level has been placed in the centre of the band-gap for
the insulating results.

Interestingly, in the Mott insulating phase of the $\rm{x}=0.25$ system Fig.
B.1 (b)-(d), there are two distinct features above the Fermi level, both of
which we associate to the upper Hubbard band. Here, the higher energy feature
is broad in nature (and is indistinguishable amongst the V-Ti hybridised
states in the main text) whereas the lower energy feature has a lower energy
shoulder which becomes suppressed with increasing $U$ until the shoulder is no
longer visible for $U=5$ eV. The (peak-to-peak) energy separation between
these two features of the upper Hubbard band remains constant as a function of
$U$. These properties of these two observed features of the upper Hubbard band
ties in well with the observed three $N=2$ (two electron) excited atomic
states calculated for bulk SrVO3 in Ref. [60]. The three $N=2$ atomic states
are: same spin singly occupied orbitals (lowest energy relative to the ground
state); different spin singly occupied orbitals; and doubly occupied orbitals
(highest energy relative to the ground state) [60]. Therefore the shoulder is
likely associated with the same spin singly occupied atomic state, the peak of
the lower feature of the upper Hubbard band relates to the different spin
singly occupied atomic state, which leaves the higher energy feature of the
upper Hubbard band being connected to the doubly occupied atomic state. In
Fig. B.1 (b) and (c) The energy difference between the two singly occupied
orbital states is much smaller than between a singly occupied orbital state
and the doubly occupied atomic state, which agrees with the atomic states seen
in bulk SrVO3 and the solution to the atomic problem of the interaction
Hamiltonian presented in Ref. [60]. All of this indicates that these spectral
function features above the Fermi level are likely associated to the $N=2$
atomic states. Analytic continuation is likely smearing out the finer details
of these atomic states in the upper Hubbard band which is why the shoulder is
being suppressed assuming the spectral weight around the $N=2$ same spin
singly occupied atomic state changes for higher $U$ values as seen in Ref.
[60].

For the $\rm{x}=0.5$ system in Fig. B.1, the upper Hubbard band is a more
broad feature similar to bulk SrVO3. The atomic states are less
distinguishable but there is a notable smoothed kink in the spectral function
above the peak upper Hubbard band energy (e.g., at about 1.8 eV in Fig. B.1
(b)) which may also indicate the $N=2$ atomic states, where the $N=2$ singly
occupied atomic states are still likely at lower energies than the doubly
occupied atomic state.

##  Appendix C External-Potential-driven Metal-Insulator Transition

Figure C.1: The effect of varying applying an effective external (gating)
potential Vg on the $\rm{x}=0.25$ and 0.5 Wannier V t2g spectral functions
(calculated from the Wannier projectors without the Ti-V hybridised states)
which pushes the insulating state through the metal-insulator transition to a
clear metallic state. The effect of Vg (0.34 V for $\rm{x}=0.25$ and 0.26 V
for $\rm{x}=0.5$) is equivalent to an excess electron doping of 0.25 per V
site (e-/site) for $\rm{x}=0.25$ and 0.025 per V site (e-/site) for
$\rm{x}=0.5$. The Fermi level has been placed in the centre of the band-gap
for the insulating results.

In the main text, we discuss the influence of applying an external potential
to the $\rm{x}=0.5$ system (of the unit cell configuration in the main text)
which pushes it from the insulating phase to the metallic one. This is shown
in Fig. C.1 along with the larger potential needed to push the $\rm{x}=0.25$
insulating system (of the unit cell configuration in the main text) to become
metallic. These DFT+DMFT results are from the calculations using $U=4.0$ eV.
These results show that the MIT from an external potential is not necessarily
confined to be around the critical doping $\rm{x}_{\rm c}$, possibly giving
more flexibility for use in Mottronic devices.

## References

  * Ha and Ramanathan [2011] S. D. Ha and S. Ramanathan, Adaptive oxide electronics: A review, Journal of Applied Physics 110, 071101 (2011).
  * Mannhart and Schlom [2010] J. Mannhart and D. G. Schlom, Oxide interfaces — an opportunity for electronics, Science 327, 1607 (2010).
  * Brahlek  _et al._ [2017] M. Brahlek, L. Zhang, J. Lapano, H.-T. Zhang, R. Engel-Herbert, N. Shukla, S. Datta, H. Paik, and D. G. Schlom, Opportunities in vanadium-based strongly correlated electron systems, MRS Communications 7, 27–52 (2017).
  * Imada  _et al._ [1998] M. Imada, A. Fujimori, and Y. Tokura, Metal-insulator transitions, Rev. Mod. Phys. 70, 1039 (1998).
  * Newns  _et al._ [1998] D. M. Newns, J. A. Misewich, C. C. Tsuei, A. Gupta, B. A. Scott, and A. Schrott, Mott transition field effect transistor, Appl. Phys. Lett. 73, 780 (1998).
  * Nekrasov  _et al._ [2005] I. A. Nekrasov, G. Keller, D. E. Kondakov, A. V. Kozhevnikov, T. Pruschke, K. Held, D. Vollhardt, and V. I. Anisimov, Comparative study of correlation effects in CaVO3 and SrVO3, Phys. Rev. B 72, 155106 (2005).
  * Zhang  _et al._ [2016] L. Zhang, Y. Zhou, L. Guo, W. Zhao, A. Barnes, H.-T. Zhang, C. Eaton, Y. Zheng, M. Brahlek, H. F. Haneef, N. J. Podraza, M. H. W. Chan, V. Gopalan, K. M. Rabe, and R. Engel-Herbert, Correlated metals as transparent conductors, Nature Materials 15, 204 (2016).
  * Mirjolet  _et al._ [2019] M. Mirjolet, F. Sánchez, and J. Fontcuberta, High carrier mobility, electrical conductivity, and optical transmittance in epitaxial srvo3 thin films, Advanced Functional Materials 29, 1808432 (2019).
  * Yoshida  _et al._ [2010] T. Yoshida, M. Hashimoto, T. Takizawa, A. Fujimori, M. Kubota, K. Ono, and H. Eisaki, Mass renormalization in the bandwidth-controlled mott-hubbard systems SrVO3 and CaVO3 studied by angle-resolved photoemission spectroscopy, Phys. Rev. B 82, 085119 (2010).
  * Aizaki  _et al._ [2012] S. Aizaki, T. Yoshida, K. Yoshimatsu, M. Takizawa, M. Minohara, S. Ideta, A. Fujimori, K. Gupta, P. Mahadevan, K. Horiba, H. Kumigashira, and M. Oshima, Self-energy on the low- to high-energy electronic structure of correlated metal SrVO3, Phys. Rev. Lett. 109, 056401 (2012).
  * Sekiyama  _et al._ [2004] A. Sekiyama, H. Fujiwara, S. Imada, S. Suga, H. Eisaki, S. I. Uchida, K. Takegahara, H. Harima, Y. Saitoh, I. A. Nekrasov, G. Keller, D. E. Kondakov, A. V. Kozhevnikov, T. Pruschke, K. Held, D. Vollhardt, and V. I. Anisimov, Mutual experimental and theoretical validation of bulk photoemission spectra of Sr1-xCaxVO3, Phys. Rev. Lett. 93, 156402 (2004).
  * Pen  _et al._ [1999] H. F. Pen, M. Abbate, A. Fujimori, Y. Tokura, H. Eisaki, S. Uchida, and G. A. Sawatzky, Electronic structure of Y1-xCaxVO3 studied by high-energy spectroscopies, Phys. Rev. B 59, 7422 (1999).
  * Laverock  _et al._ [2013] J. Laverock, B. Chen, K. E. Smith, R. P. Singh, G. Balakrishnan, M. Gu, J. W. Lu, S. A. Wolf, R. M. Qiao, W. Yang, and J. Adell, Resonant soft-x-ray emission as a bulk probe of correlated electron behavior in metallic SrxCa1-xVO3, Phys. Rev. Lett. 111, 047402 (2013).
  * Nekrasov  _et al._ [2006] I. A. Nekrasov, K. Held, G. Keller, D. E. Kondakov, T. Pruschke, M. Kollar, O. K. Andersen, V. I. Anisimov, and D. Vollhardt, Momentum-resolved spectral functions of SrVO3 calculated by LDA+DMFT, Phys. Rev. B 73, 155112 (2006).
  * Kotliar  _et al._ [2006] G. Kotliar, S. Y. Savrasov, K. Haule, V. S. Oudovenko, O. Parcollet, and C. A. Marianetti, Electronic structure calculations with dynamical mean-field theory, Rev. Mod. Phys. 78, 865 (2006).
  * Byczuk  _et al._ [2007] K. Byczuk, M. Kollar, K. Held, Y.-F. Yang, I. A. Nekrasov, T. Pruschke, and D. Vollhardt, Kinks in the dispersion of strongly correlated electrons, Nat. Phys. 3, 168 (2007).
  * Tomczak  _et al._ [2012] J. M. Tomczak, M. Casula, T. Miyake, F. Aryasetiawan, and S. Biermann, Combined GW and dynamical mean-field theory: Dynamical screening effects in transition metal oxides, Europhys. Lett. 100, 67001 (2012).
  * Georges  _et al._ [1996] A. Georges, G. Kotliar, W. Krauth, and M. J. Rozenberg, Dynamical mean-field theory of strongly correlated fermion systems and the limit of infinite dimensions, Rev. Mod. Phys. 68, 13 (1996).
  * Boehnke  _et al._ [2016] L. Boehnke, F. Nilsson, F. Aryasetiawan, and P. Werner, When strong correlations become weak: Consistent merging of $GW$ and DMFT, Phys. Rev. B 94, 201106 (2016).
  * Zhong  _et al._ [2015] Z. Zhong, M. Wallerberger, J. M. Tomczak, C. Taranto, N. Parragh, A. Toschi, G. Sangiovanni, and K. Held, Electronics with Correlated Oxides: SrVO3/SrTiO3 as a Mott Transistor, Phys. Rev. Lett. 114, 246401 (2015).
  * Yoshimatsu  _et al._ [2010] K. Yoshimatsu, T. Okabe, H. Kumigashira, S. Okamoto, S. Aizaki, A. Fujimori, and M. Oshima, Dimensional-crossover-driven metal-insulator transition in SrVO3 ultrathin films, Phys. Rev. Lett. 104, 147601 (2010).
  * Gu  _et al._ [2014] M. Gu, S. A. Wolf, and J. Lu, Two-dimensional Mott insulators in SrVO3 ultrathin films, Adv. Mater. Interfaces 1, 1300126 (2014).
  * Bhandary  _et al._ [2016] S. Bhandary, E. Assmann, M. Aichhorn, and K. Held, Charge self-consistency in density functional theory combined with dynamical mean field theory: $k$-space reoccupation and orbital order, Phys. Rev. B 94, 155131 (2016).
  * Schüler  _et al._ [2018] M. Schüler, O. E. Peil, G. J. Kraberger, R. Pordzik, M. Marsman, G. Kresse, T. O. Wehling, and M. Aichhorn, Charge self-consistent many-body corrections using optimized projected localized orbitals, J. Phys.: Conden. Matter 30, 475901 (2018).
  * Hampel  _et al._ [2020] A. Hampel, S. Beck, and C. Ederer, Effect of charge self-consistency in $\mathrm{DFT}+\mathrm{DMFT}$ calculations for complex transition metal oxides, Phys. Rev. Research 2, 033088 (2020).
  * Beck  _et al._ [2018] S. Beck, G. Sclauzero, U. Chopra, and C. Ederer, Metal-insulator transition in CaVO3 thin films: Interplay between epitaxial strain, dimensional confinement, and surface effects, Phys. Rev. B 97, 075107 (2018).
  * Sclauzero  _et al._ [2016] G. Sclauzero, K. Dymkowski, and C. Ederer, Tuning the metal-insulator transition in $d^{1}$ and $d^{2}$ perovskites by epitaxial strain: A first-principles-based study, Phys. Rev. B 94, 245109 (2016).
  * James  _et al._ [2021a] A. D. N. James, M. Aichhorn, and J. Laverock, Quantum confinement induced metal-insulator transition in strongly correlated quantum wells of ${\mathrm{SrVO}}_{3}$ superlattices, Phys. Rev. Research 3, 023149 (2021a).
  * Laverock  _et al._ [2017] J. Laverock, M. Gu, V. Jovic, J. W. Lu, S. Wolf, R. M. Qiao, W. Yang, and K. Smith, Nano-engineering of electron correlation in oxide superlattices, Nano Futures 1, 031001 (2017).
  * Gu  _et al._ [2018] M. Gu, S. A. Wolf, and J. Lu, Transport phenomena in SrVO3/SrTiO3 superlattices, J. Phys. D: Appl. Phys. 51, 10LT01 (2018).
  * Kobayashi  _et al._ [2015] M. Kobayashi, K. Yoshimatsu, E. Sakai, M. Kitamura, K. Horiba, A. Fujimori, and H. Kumigashira, Origin of the anomalous mass renormalization in metallic quantum well states of strongly correlated oxide SrVO3, Phys. Rev. Lett. 115, 076801 (2015).
  * Kobayashi  _et al._ [2017] M. Kobayashi, K. Yoshimatsu, T. Mitsuhashi, M. Kitamura, E. Sakai, R. Yukawa, M. Minohara, A. Fujimori, K. Horiba, and H. Kumigashira, Emergence of quantum critical behavior in metallic quantum-well states of strongly correlated oxides, Sci. Rep. 7, 16621 (2017).
  * Rossella  _et al._ [2007] F. Rossella, P. Galinetto, G. Samoggia, V. Trepakov, and L. Jastrabik, Photoconductivity and the structural phase transition in srtio3, Solid State Communications 141, 95 (2007).
  * Li  _et al._ [2019] R. Li, C. Zhang, J. Liu, J. Zhou, and L. Xu, A review on the electrical properties of doped srtio3 as anode materials for solid oxide fuel cells, Materials Research Express 6, 102006 (2019).
  * Xu  _et al._ [2022] Y. Xu, Y. Liang, Q. He, R. Xu, D. Chen, X. Xu, and H. Hu, Review of doping srtio3 for photocatalytic applications, Bulletin of Materials Science 46, 6 (2022).
  * Kleemann  _et al._ [2020] W. Kleemann, J. Dec, A. Tkach, and P. M. Vilarinho, SrTiO3 \- Glimpses of an Inexhaustible Source of Novel Solid State Phenomena, Condensed Matter 5, 10.3390/condmat5040058 (2020).
  * Shi  _et al._ [2020] X.-L. Shi, H. Wu, Q. Liu, W. Zhou, S. Lu, Z. Shao, M. Dargusch, and Z.-G. Chen, SrTiO3-based thermoelectrics: Progress and challenges, Nano Energy 78, 105195 (2020).
  * Yan  _et al._ [2018] H. Yan, Z. Zhang, S. Wang, and K. Jin, Review of photoresponsive properties at SrTiO3-based heterointerfaces, Chinese Physics B 27, 117804 (2018).
  * Pai  _et al._ [2018] Y.-Y. Pai, A. Tylan-Tyler, P. Irvin, and J. Levy, Physics of SrTiO3-based heterostructures and nanostructures: a review, Reports on Progress in Physics 81, 036503 (2018).
  * Christensen  _et al._ [2019] D. V. Christensen, F. Trier, W. Niu, Y. Gan, Y. Zhang, T. S. Jespersen, Y. Chen, and N. Pryds, Stimulating Oxide Heterostructures: A Review on Controlling SrTiO3-Based Heterointerfaces with External Stimuli, Advanced Materials Interfaces 6, 1900772 (2019), https://onlinelibrary.wiley.com/doi/pdf/10.1002/admi.201900772 .
  * Miyasaka  _et al._ [2000] S. Miyasaka, T. Okuda, and Y. Tokura, Critical Behavior of Metal-Insulator Transition in ${\mathrm{La}}_{1-\mathit{x}}{\mathrm{Sr}}_{\mathit{x}}{\mathrm{VO}}_{3}$, Phys. Rev. Lett. 85, 5388 (2000).
  * Takahashi  _et al._ [2022] K. S. Takahashi, Y. Tokura, and M. Kawasaki, Metal-insulator transitions in dimensionality controlled ${\mathrm{La}}_{1-x}{\mathrm{Sr}}_{x}{\mathrm{VO}}_{3}$ films, APL Materials 10, 111114 (2022).
  * Hong  _et al._ [2002] K. Hong, S.-H. Kim, Y.-J. Heo, and Y.-U. Kwon, Metal-insulator transitions of SrTi1-xVxO3 solid solution system, Solid State Communications 123, 305 (2002).
  * Kanda  _et al._ [2021] T. Kanda, D. Shiga, R. Yukawa, N. Hasegawa, D. K. Nguyen, X. Cheng, R. Tokunaga, M. Kitamura, K. Horiba, K. Yoshimatsu, and H. Kumigashira, Electronic structure of $\mathrm{Sr}{\mathrm{Ti}}_{1\text{$-$}x}{\mathrm{V}}_{x}{\mathrm{O}}_{3}$ films studied by in situ photoemission spectroscopy: Screening for a transparent electrode material, Phys. Rev. B 104, 115121 (2021).
  * Gu  _et al._ [2013] M. Gu, S. A. Wolf, and J. Lu, Metal-insulator transition in SrTi1-xVxO3 thin films, Applied Physics Letters 103, 223110 (2013).
  * Tsuiki  _et al._ [1983] H. Tsuiki, K. Kitazawa, and K. Fueki, The Donor Level of V4+ and the Metal-Nonmetal Transition in SrTi1-xVxO3, Japanese Journal of Applied Physics 22, 590 (1983).
  * Mantry and Sarun [2019] S. P. Mantry and P. Sarun, Investigation on the impedance spectroscopy and electrical conduction mechanism in SrTi1-xVxO3 ($x=0.00\leq x\leq 0.15$) ceramics, Journal of Materials Science: Materials in Electronics 30, 6795 (2019).
  * Itaka  _et al._ [2005] K. Itaka, T. Wakisaka, T. Mihara, S. Yaginuma, Y. Matsumoto, and H. Koinuma, Sharp metal-insulator transition in SrTi1-xVxO3 thin films on SrTiO3 substrates, Thin Solid Films 486, 222 (2005), wOE-11, 2004.
  * Sturge [1968] M. Sturge, The Jahn-Teller Effect in Solids (Academic Press, 1968) pp. 91–211.
  * [50] J. K. Dewhurst, S. Sharma, L. Nordström, F. Cricchio, O. Granas, and E. K. U. Gross, http://elk.sourceforge.net/.
  * Perdew  _et al._ [1996] J. P. Perdew, K. Burke, and M. Ernzerhof, Generalized gradient approximation made simple, Phys. Rev. Lett. 77, 3865 (1996).
  * Parcollet  _et al._ [2015] O. Parcollet, M. Ferrero, T. Ayral, H. Hafermann, I. Krivenko, L. Messio, and P. Seth, TRIQS: A toolbox for research on interacting quantum systems, Comp. Phys. Commun. 196, 398 (2015).
  * James  _et al._ [2021b] A. D. N. James, E. I. Harris-Lee, A. Hampel, M. Aichhorn, and S. B. Dugdale, Wave functions, electronic localization, and bonding properties for correlated materials beyond the Kohn-Sham formalism, Phys. Rev. B 103, 035106 (2021b).
  * Aichhorn  _et al._ [2016] M. Aichhorn, L. Pourovskii, P. Seth, V. Vildosola, M. Zingl, O. E. Peil, X. Deng, J. Mravlje, G. J. Kraberger, C. Martins, M. Ferrero, and O. Parcollet, TRIQS/DFTTools: A TRIQS application for ab initio calculations of correlated materials, Comp. Phys. Commun. 204, 200 (2016).
  * Seth  _et al._ [2016] P. Seth, I. Krivenko, M. Ferrero, and O. Parcollet, TRIQS/CTHYB: A continuous-time quantum Monte Carlo hybridisation expansion solver for quantum impurity problems, Comp. Phys. Commun. 200, 274 (2016).
  * Aichhorn  _et al._ [2009] M. Aichhorn, L. Pourovskii, V. Vildosola, M. Ferrero, O. Parcollet, T. Miyake, A. Georges, and S. Biermann, Dynamical mean-field theory within an augmented plane-wave framework: Assessing electronic correlations in the iron pnictide LaFeAsO, Phys. Rev. B 80, 085101 (2009).
  * Kraberger  _et al._ [2017] G. J. Kraberger, R. Triebl, M. Zingl, and M. Aichhorn, Maximum entropy formalism for the analytic continuation of matrix-valued Green’s functions, Phys. Rev. B 96, 155128 (2017).
  * Georgescu [2023] A. B. Georgescu, Cu-doped Pb10(PO4)6O, and V doped SrTiO3 – a tutorial on electron-crystal lattice coupling in insulating materials with transition metal dopants (2023), arXiv:2308.07295 .
  * Laverock  _et al._ [2015] J. Laverock, J. Kuyyalil, B. Chen, R. P. Singh, B. Karlin, J. C. Woicik, G. Balakrishnan, and K. E. Smith, Enhanced electron correlations at the SrxCa1-xVO3 surface, Phys. Rev. B 91, 165123 (2015).
  * Bauernfeind  _et al._ [2017] D. Bauernfeind, M. Zingl, R. Triebl, M. Aichhorn, and H. G. Evertz, Fork tensor-product states: Efficient multiorbital real-time dmft solver, Phys. Rev. X 7, 031013 (2017).
  * Nguyen  _et al._ [2022] T. H. Y. Nguyen, D. A. Le, and A. T. Hoang, Anderson localization in the Anderson-Hubbard model with site-dependent interactions, New Journal of Physics 24, 053054 (2022).
  * Weh  _et al._ [2021] A. Weh, Y. Zhang, A. Östlin, H. Terletska, D. Bauernfeind, K.-M. Tam, H. G. Evertz, K. Byczuk, D. Vollhardt, and L. Chioncel, Dynamical mean-field theory of the anderson-hubbard model with local and nonlocal disorder in tensor formulation, Phys. Rev. B 104, 045127 (2021).


# Decentralized Collaborative Pricing and Shunting for Multiple EV Charging
Stations Based on Multi-Agent Reinforcement Learning

1st Tianhao Bu  Glory Engineering & Tech Co., LTD
Shanghai, China
<EMAIL_ADDRESS>2nd Hang Li⋆ School of Electronic
Information
and Electrical Engineering
Shanghai Jiao Tong University
Shanghai,China
<EMAIL_ADDRESS>
⋆Corresponding author  3rd Guojie Li  School of Electronic Information
and Electrical Engineering
Shanghai Jiao Tong University
Shanghai,China
<EMAIL_ADDRESS>
###### Abstract

The extraordinary electric vehicle (EV) popularization in the recent years has
facilitated research studies in alleviating EV energy charging demand.Previous
studies primarily focused on the optimizations over charging stations’ (CS)
profit and EV users’ cost savings through charge/discharge scheduling events.
In this work, the random behaviours of EVs are considered, with EV users’
preferences over multi-CS characteristics modelled to imitate the potential CS
selection disequilibrium. A price scheduling strategy under decentralized
collaborative framework is proposed to achieve EV shunting in a multi-CS
environment, while minimizing the charging cost through multi-agent
reinforcement learning. The proposed problem is formulated as a Markov
Decision Process (MDP) with uncertain transition probability.

###### Index Terms:

electric vehicle price scheduling, electric vehicle shunting, multi-agent
reinforcement learning, centralised training decentralised execution, markov
decision process

##  I INTRODUCTION

Figure 1: Decentralized collaborative framework for a single CS.

According to the report from the House of Common Library [1], UK has proposed
and adopted two sets of energy consumption strategies to reach the net zero
policy by 2050. The vehicle carbon emissions have occupied the largest
proportion of UK emission of 23‰ in 2022. Within the same year, UK achieved
the second largest electric vehicle (EV) sales in Europe [2]. The convincing
results further propelled the policies over petrol vehicle prohibition and EV
mandate, in turn facilitated the EV popularization thus the ever-growing
charging demand. Correspondingly, developing solutions to alleviate the demand
crisis and reduce charging peak load has gained noticeable attentions in the
recent years.
The nature of the aforementioned problems requires sequential decision-making
approaches, thus implying them into Markov Decision Process (MDP) [3] is
favoured. Reinforcement learning (RL) as one of the main stream machine
learning paradigms has been a hot topic within this domain, due to its
capability of learning optimal policies in complicated environments in a
model-free manner. RL is frequently applied in managing the EV charging,
discharging schedulings to optimize charging costs [4] , [5]; In [5], RL
algorithm is used to learn device based MDP models, with the implementation of
a Q-table to estimate the Q function. At the other hand, [4] uses recurrent
deep deterministic policy gradient (RDDPG) algorithm, which proven to have
great scalability in solving large-scaled MDP problems.
While the above methods focus on single agent learning problems, the real-
world environment involves multi-EVs charging under single or multi-CSs. The
scenarios where multi-CS are participated in, the inherent competitive and
cooperative natures among CSs should also be considered. To adapt with the
corresponding increased complexity, the emergence of multi-agent reinforcement
learning (MARL) based approaches opens up a new realm worth investigation [6],
[7], [8], [9]. MARL based methods can be generalized into two categories 1)
centralized execution methods [8], [9], 2) decentralized execution methods
[6], [7]. For centralized training and decentralized execution (CTDE), [6]
proposes an actor-attention-critic approach to boost reward optimization
performance through an attention critic allowing dynamic agents selections at
each training time point; and in [7], a decentralized collaborative framework
is proposed to account situations where CS has various types of charging
piles, i.e. normal charging (NC) and fast charging (FC).
Although the aforementioned methods are sufficient in providing charging cost
optimization improvements, there weren’t many works on modelling and, or
guiding EV users’ characteristics. One of the branches conducts research
within the direction of EV charging navigation [10], [11], where the primary
objective is to suggest users with the shortest routes between EVs and CSs,
while achieving simultaneous minimization over EV users’ traveling times and
charging costs. In [10], a deterministic shortest charging route model (DSCRM)
is utilized to extract state features from the stochastic traffic, and RL is
applied to obtain an optimal strategy for the navigation problem with a
practical case study in Xi’an’s real-sized zone. Whereas in [11], the features
are extracted through graph convolutional transformation and learnt through
graph reinforcement learning (GRL), the optimization goal is achieved with the
aid of a modified rainbow algorithm. Other directions consider modelling the
EV users’ preferences [12], [13]. [12] implements a user preference and price
based charging algorithm, experimented in the UCLA parking lots; and [13]
proposed a modified probability model.
In this work, the decentralized collaborative framework is extended to include
multi-CSs, the CS characteristics are represented in an asymmetric manner,
i.e. diversified charging prices at different time periods, varied charging
pile sizes. To study the EV users’ behaviours, a probabilistic EV user
preference model is built accounting real time charging prices, EV to CS
distances and CS sizes (pile numbers) . Hence, a normalized linear distance
model is also developed along with the random arrivals of the EVs. Finally,
since the probabilistic EV user preferences can potentially result in an
asymmetric selection over a certain CS, a price scheduling strategy is
proposed to allow CSs to compete over prices and attract EV users in real time
through shared CS information, and eventually achieve EV shunting, relieving
the stress of the potential space congestion and charging demands over certain
individual CS.

The main contributions of this paper are summarized as follows:

  1. 1.

A multi-CS based decentralized collaborative framework is proposed, where each
pile with an EV attached is treated as an agent, which is able to execute
local and independent charging, discharging actions; the charging problem is
formulated as MDP. A Value-Decomposition Network (VDN) [14] MARL algorithm is
utilized to perform charging cost optimization.

  2. 2.

A probabilistic EV user preference model is designed, which considers real
time prices, EV-CS distances and CS sizes to predict EV users selections over
various CSs; each term is assigned with a tuned weight, together with a linear
distance model developed to normalize and quantify the distance parameter.

  3. 3.

A real time dynamic price scheduling strategy is proposed, by utilizing the
shared information between the stations, individual CS is able to influence EV
users’ selections through live price competition; the reduced price for the
influenced EV is reclaimed in the remaining time steps before the departure to
ensure CS profit, while achieving EV shunting.

  4. 4.

Through VDN and Q-mix [15] performance evaluations, the results show that our
price scheduling strategy is sufficient in producing substantial average
reward improvements and faster episode convergence over the baseline.

The rest of the paper is structured as follows. section II presents the
concept of single CS based decentralized collaborative framework. section III
proposes the multi-CS based collaborative framework, followed with the
modified probabilistic EV user preference model definition, and the MDP
problem formulation. section IV explains the proposed price scheduling
strategy procedure and the training phase. In section V, rigorous evaluations
were carried out to validate the effectiveness of the proposed method,
together with the corresponding optimisation performances. Finally, section VI
concludes our work.

##  II BACKGROUND

###  II-A Decentralized collaborative framework (single CS)

Figure 2:  Multi(dual)-CS based decentralized collaborative framework

Under a single CS based decentralized collaborative framework, the CS contains
multiple charging piles, including both NC and FC, as shown in Figure 1. Once
EV arrives and attaches itself to an allocated pile, its information, i.e.
current state of charge (SOC), time of arrival $t_{a}$, estimate time of
departure $t_{d}$ and the charging demand $e_{d}$ will be transmitted via the
pile and CS to the data operator; presented as data signals.Since each pile
has the ability to make local and independent charging, discharging decisions
(as forms of energy flow exchanges), data signals also embed these decisions
together with the CS real time charging prices. In the single CS scenario, the
data operator’s task is to gather these data and send control signals to
ensure the accumulated independent charging decisions will not cause CS
capacity overload , this is implemented through live monitoring of the total
CS power demand at each time step, in the same manner as the procedure
proposed in [4]. In our work, the single CS based framework is extended to a
dual-CS decentralized collaborative framework, as the result, the data
operator also acts as a medium between the CSs to exchange information of the
charging prices, real time charged EV number of each CS, allowing the
proceeding of the proposed price ”competition” algorithm; the detailed
descriptions are presented in section IV.

##  III ENVIRONMENT FORMULATION

###  III-A Multi(dual)-CS based collaborative framework

From the single-CS based framework mentioned in section II, we propose a
multi-CS based collaborative framework by extending the single station into
dual stations, as depicted in Figure 2. Under the new framework, CS 0 and CS 1
are differentiated by their 1) Charging pile number $N_{z}$ , with
$z\in[0..1]$. 2) Real time charging price $p_{t,z}$. The data operator now
carries additional responsibility of integrating and exchanging the pile
number and price information between the CSs. The proposed framework inherits
the same concepts as its baseline, in which each arrived EV, equivalent as a
pile under operation, is treated as an agent, moreover, the optimization
objective is to minimize the accumulated charging cost $C_{acc}$, through each
individual pile’s charging decisions, expressed as

$C_{acc}=\sum_{z}\sum_{t}p_{t,z}$ (1)

###  III-B Probabilistic EV user preference model

Inspired by the works from [16] and [13], we propose a modified probabilistic
model to imitate EV users’ selection preferences over the CSs. The following
parameters being used are adopted from [16].
The selection probability $Pr_{i,z,t_{a}^{i}}$ of the $i^{th}$ arriving EV
selecting CS $z$, at its arrival time $t_{a}^{i}$ can be formally defined as
follows:

$Pr_{i,z,t_{a}^{i}}=\frac{U_{i,t_{a}^{i},z}}{\sum_{z=0}^{1}U_{i,t_{a}^{i},z}}$
(2)

where $U_{i,t_{a}^{i},z}$ is the EV user utility function, in turn is directly
proportional to the attractiveness $Att_{i,t_{a}^{i},z}$ of CS $z$ on the
$i^{th}$ EV user . The user utility function can be expressed as

$U_{i,t_{a}^{i},z}=\frac{Att_{i,t_{a}^{i},z}}{T_{i,z}}$ (3)

where $T_{i,z}$ defines the closeness between EV $i$ and CS $z$, represented
as

$T_{i,z}=t_{c}+t_{i,z}$ (4)

with $t_{c}$ being the average EV charging time and $t_{i,z}$ being the
travelling time from EV $i$ to CS $z$.
In a practical scenario, we consider four major factors that have significant
influences on EV users’ decisions: 1) Charging pile number $N_{z}$. By
considering EV users’ risk stop loss characteristics, the arrived users’
selections would lean towards the CS with higher pile number, due to the
likelihood of having more available spaces, and the smaller chances of
discovering fully occupied station at arrival, which results in longer
distance travelled, towards inferior charging locations. 2) Arrival time
charging price $p_{t_{a}^{i},z}$. EV users generally would prefer cheaper
prices, hence, the CS that offers cheaper prices will gain more attractions
from the users. 3) EV $i$ and CS $z$ travel distance $d_{i,z}$. Longer travel
distance will make CS less attractive to the EV users. 4) EV arrival state of
charge $SOC_{a}$. High $SOC_{a}$ will increase CS $z$’s attractiveness.
Furthermore, we define a SOC threshold $SOC_{th}$, when $SOC_{a}$ is less than
$SOC_{th}$, the EV is urgent to get charged, leading to less attraction under
large $d_{i,z}$, and vice versa. This term will be discussed in more details
in the upcoming session.
Based on the aforementioned factors, the attractiveness $Att_{i,t_{a}^{i},z}$
can be formulated as a linear weighing function:

$Att_{i,t_{a}^{i},z}=w_{0}\cdot N_{z}-w_{1}\cdot p_{z,t_{a}^{i}}-w_{2}\cdot
d_{i,z}+w_{3}\cdot(SOC_{a}-SOC_{th})$ (5)

$w_{j}$ with $j\in[0..3]$ are the weighing coefficients, the negative sign
indicates degrading effect on the user attraction.

Figure 3: Linear EV-CS distance model

###  III-C Distance normalization model

To illustrate the concept of utilizing distance normalization, we first
propose a linear EV-CS distance model, as shown in Figure 3. $d_{i,0}$ and
$d_{i,1}$ are the corresponding distances between EV $i$ and CS 0, 1, and they
are assumed to be the shortest routes by default. Recalling from equation (5),
high $SOC_{a}$ would result in higher station attraction to EV users, and the
SOC is related to the EV-CS distances. But due to the unit differences (SOC
has values ranging from 0 to 1 while the actual distance $d_{i,z}$ is
significantly larger in comparison), they are separated and assigned with
different weighing coefficients $w_{2}$ and $w_{3}$, in turn increases the
attraction function complexity and obscures the relationship between them.
Hence, we divide the original $d_{i,z}$ by the total distance from EV $i$ to
both stations, and obtain the normalized distance $\widetilde{d_{i,z}}$,
expressed as

$\widetilde{d_{i,z}}=\frac{d_{i,z}}{d_{i,0}+d_{i,1}},z\in[0..1]$ (6)

As a result, the original attractiveness $Att_{i,t_{a}^{i},z}$ can be modified
into $\widetilde{Att}_{i,t_{a}^{i},z}$, as shown the following formula:

$\widetilde{Att}_{i,t_{a}^{i},z}=w_{0}\cdot N_{z}-w_{1}\cdot
p_{t_{a}^{i},z}-w_{2}\cdot(\widetilde{d_{i,z}}+SOC_{th}-SOC_{a})$ (7)

From equation (26), $w_{3}$ is eliminated as $\widetilde{d_{i,z}}$, $SOC_{a}$
and $SOC_{th}$ now share a highly overlapped value range. When the arrival SOC
exceeds the threshold, the observation indicates the EV having sufficient
energy to travel relative long distances, in turn alleviates the attraction
degradation caused from $\widetilde{d_{i,z}}$, and vice versa.
Additionally, the traveling time $t_{i,z}$ can be calculated using
$\widetilde{d_{i,z}}$:

$t_{i,z}=\frac{\widetilde{d_{i,z}}}{\widetilde{v_{i}}}$ (8)

where $\widetilde{v_{i}}$ is the normalized vehicle velocity of EV $i$.

###  III-D Markov decision process

As a discrete-time stochastic control process, MDP is sufficient in modeling
our charging,discharging decision problems, due to the sequential feature of
the decision making process. MDP consists of four components: 1) State. It
describes the agent’s (arrived EV, or the pile under operation) current
situation, and the future ’next’ state only depends on the present state
,owing to the MDP’s definition. 2) Action. The charging, discharging decisions
made by the NC, FC piles, under the current state time step. 3) State
transition function. The function that evaluates the probability in which the
action under current state causing the future ’next’ state transition. 4)
Reward function. The optimization goal when the action causes current state
transition to the future ’next’ state, and it’s directly related to the
accumulated charging cost.
In this section, we use MDP to model our decision problem, with detailed
definitions presented in the following subsections.

####  III-D1 State

The state at time step $t$ can be defined as
$s_{i,t,z}=(s_{i,t,z}^{pile},s_{t,z}^{station})$, where $s_{i,t,z}^{pile}$ is
the state of pile $i$ in CS $z$ at time $t$, which can be formally defined as

$s_{i,t,z}^{pile}=(SOC_{i,t,z},P^{max}_{i,t,z},P^{min}_{i,t,z},t^{stay}_{i,t,z})$
(9)

where $t^{stay}_{i,t,z}$ is the remaining staying time of the $i^{th}$ EV in
CS $z$ at time step $t$, defined as

$t^{stay}_{i,t,z}=t^{d}_{i,t,z}-t^{cur}$ (10)

with $t^{d}_{i,t,z}$ being the time of departure, and $t^{cur}$being the
current time step.
Furthermore, $s_{t,z}^{station}$ is the state of CS $z$ at time $t$, and can
be expressed as

$s_{t,z}^{station}=(t^{cur},p_{t,z},eme^{total}_{t,z})$ (11)

where $eme^{total}_{t,z}$ is CS $z$’s total emergency at time $t$, can also in
turn be represented as

$eme^{total}_{t,z}=\sum_{i}^{N_{i,t,z}}eme_{i,t,z}$ (12)

with $N_{i,t,z}$ being the number of EVs attached to CS $z$ at time $t$, and
$eme_{i,t,z}$ being the charge emergence of the individual EV $i$ at CS $z$
and time step $t$, as shown in (13).

$eme_{i,t,z}=\begin{cases}\frac{(SOC_{i,t,z}^{need}-SOC_{i,t,z})\cdot
C_{i,z}}{t^{stay}_{i,t,z}}&\text{if $SOC_{i,t,z}$ $<$ $SOC_{i,t,z}^{need}$}\\\
0&\text{otherwise}\par\end{cases}$ (13)

$SOC_{i,t,z}^{need}$ is the SOC needed for EV $i$ at time $t$ charging at
station $z$ before its departure.

####  III-D2 Action

In this section we define the $i^{th}$ pile’s action in CS $z$ at time step
$t$ as $a_{i,t,z}$, within a range of $[-1,1]$. $a_{i,t,z}$ represents the
corresponding charging,discharging decision made by the pile. Thus, the real
time charging power $P^{real}_{i,t,z}$ is expressed as

$P^{real}_{i,t,z}=\frac{(a_{i,t,z}+1)}{2}*(P^{max}_{i,t,z}-P^{min}_{i,t,z})+P^{min}_{i,t,z}$
(14)

where $P^{max}_{i,t,z}$ and $P^{min}_{i,t,z}$ are the corresponding maximum
and minimum pile charging,discharging powers at CS $z$ at time $t$.As a
result, we defined the total charging power $P^{total}_{t,z}$ in CS $z$, time
step $t$ as

$P^{total}_{t,z}=\sum_{i}P^{real}_{i,t,z}$ (15)

####  III-D3 State transition

To model the probability of state transition from $s_{t}$ to $s_{t+1}$ under
the influence of $a_{t}$, we formulate the state transition as

$s_{t+1}=f(s_{t},a_{t},\omega_{t})$ (16)

$f(\cdot)$ is the state transition function, and the term $\omega_{t}$ is
defined to indicate the uncertainty of the state transition.

####  III-D4 Reward function

The station reward $\mathbf{r}^{station}_{t,z}$ for CS $z$ at time step $t$
can be formally defined as

$\mathbf{r}^{station}_{t,z}=P^{total}_{t,z}\cdot p_{t,z}$ (17)

To ensure the total charging power $P^{total}_{t,z}$ does not exceed the
maximum power limit $P^{max}_{z}$ for each CS $z$, we define a punishment term
$\mathbf{r}^{punish}_{t,z}$ as

$\mathbf{r}^{punish}_{t,z}=\begin{cases}P^{max}_{z}-P^{total}_{t,z}&\text{if
$P^{max}_{z}$ $<$ $P^{total}_{t,z}$}\\\ 0&\text{otherwise}\end{cases}$ (18)

Therefore, the total reward at each time step $t$ is expressed as

$\mathbf{r}_{t}=\sum_{z=0}^{1}\mathbf{r}^{c}_{t,z}+A\cdot\mathbf{r}^{punish}_{t,z}$
(19)

Where $A$ is a scaling factor.

1

2

input : state $s_{i,t,z}=(s_{i,t,z}^{pile},s_{t,z}^{station})$

3

output : VDN network parameters

4 Randomly initialize VDN network parameter $\theta$.

5 Copy VDN network parameter to target network $\theta^{\prime}=\theta$.

6 for  _episode $m=1:M$ _ do

7 for  _$t=t_{start}:t_{end}$_ do

8 for  _CS $z$ = CS 0 : CS 1_ do

9 for  _pile $i$ = pile 0 : pile $n$_ do

10 obtain EV $i$’s state $s_{t,z}=(s_{i,t,z}^{pile},s_{t,z}^{station})$

11 if  _$p_{i-1,t^{cur}-1,0}\neq 0$ , $p_{i-1,t_{cur}-1,1}\neq 0$_ then

12 $p_{i-1,t-1,0}+\alpha$, $p_{i-1,t-1,1}+\beta$.

13 if  _$N_{t^{cur},0}$ $<$ $N_{t^{cur},1}$_ then

14 $p_{i,t^{cur}+1,z}=p_{t^{cur},0}-\alpha$

15 if  _$N_{t^{cur},0}$ $>$ $N_{t^{cur},1}$_ then

16 $p_{i,t^{cur}+1,z}=p_{t^{cur},0}-\beta$

17 select action $a_{i,t,z}$ based on $\epsilon-$greedy policy

18 calculate reward $\mathbf{r}^{station}_{t,z}=P^{total}_{t,z}\cdot p_{t,z}$

19 sample transition $<s_{i,t,z},a_{i,t,z},r_{t,z},s_{i,t+1,z}>$ in $B_{1}$

20

21 calculate total reward
$\mathbf{r}_{t}=\sum_{z=0}^{1}\mathbf{r}^{c}_{t,z}+A\cdot\mathbf{r}^{punish}_{t,z}$

22

23sample mini batch $<s_{i,z},a_{z},r_{z},s_{i+1,z}>$ from $B_{1}$

24 set $y^{tot}=r_{i,z}+\gamma\cdot
max_{\mathbf{a}^{\prime}}Q_{tot}({\boldsymbol{\tau}}^{\prime},\mathbf{a}^{\prime},s^{\prime}_{i,z};\theta^{\prime})$

25 update VDN network by minimizing loss
$L_{\theta}=\sum_{i}\sum_{z}(y_{i,z}^{tot}-Q_{tot}({\boldsymbol{\tau}},\mathbf{a},s_{i,z};\theta))^{2}$

Algorithm 1 EV price scheduling training procedure

##  IV PROPOSED PRICE SCHEDULING STRATEGY

###  IV-A Price scheduling

Figure 4: EV distributions over price scheduling and baseline methods, under
each scenario

In this section, we propose the novel price scheduling strategy to achieve EV
shunting through multi-CS price ”competition” at each simultaneous time step.
At $t^{cur}$, the data operator receives the pile occupation information
$N_{t^{cur},0}$ and $N_{t^{cur},1}$ from each CS, if the pile occupation is
imbalanced,$N_{t^{cur},0}\neq N_{t^{cur},1}$; the data operator will command
the CS with lower attached EV to reduce the charging price by a sufficient
amount to attract EV at step $t^{cur}+1$, namely $p_{i,t^{cur}+1,z}$ , as
shown in equation 20.

$p_{i,t^{cur}+1,z}=\begin{cases}p_{t^{cur},0}-\alpha_{t^{cur}}&\text{if
$N_{t^{cur},0}$ $<$ $N_{t^{cur},1}$}\\\ p_{t^{cur},1}-\beta_{t^{cur}}&\text{if
$N_{t^{cur},0}$ $>$ $N_{t^{cur},1}$}\\\
p_{t^{cur},z}&\text{otherwise}\par\end{cases}$ (20)

The terms $\alpha_{t^{cur}}$ and $\beta_{t^{cur}}$ are the corresponding price
reductions at each CS and the current time $t^{cur}$.
Moreover, $p_{i,t,z}$ is a subset of CS charging price $p_{t,z}$, it is only
implied to those attracted EVs at time step $t$, and will be charged back at
the next step $t+1$, in turn can be expressed as

$p_{i,t+1,z}=\begin{cases}p_{t-1,0}+\alpha_{t-1}&\text{if $z=0$,
$p_{i-N_{t-1,0}^{att},t-1,0}\neq 0$}\\\ p_{t-1,1}+\beta_{t-1}&\text{if
$z=1$,$p_{i-N_{t-1,1}^{att},t-1,1}\neq 0$}\par\end{cases}$ (21)

where $N_{t-1,0}^{att}$ and $N_{t-1,1}^{att}$ are the numbers of price
attracted EVs to CS 0 and 1 at time step $t-1$.

###  IV-B Training phase

In our work, we use VDN [14] MARL algorithm,a Q-learning based approach to
resolve the proposed multi-agent reinforcement learning problem. The primary
task of the network is to learn the total action value function
$Q_{tot}({\boldsymbol{\tau}},\mathbf{a})$; where
${\boldsymbol{\tau}}=(\tau_{1,0},...,\tau_{n,1})$ is observation history, with
$\tau_{i,z}$ being the joint action observation history of agent $i$ from CS
$z$ and $\mathbf{a}=(a_{1,0},...,a_{n,1})$ is the action set.
$Q_{tot}({\boldsymbol{\tau}},\mathbf{a})$ can be represented as a summation of
each value function $Q_{i,z}(\tau_{i,z},a_{i,t,z})$:

$Q^{tot}({\boldsymbol{\tau}},\mathbf{a})=\sum_{i}Q_{i,z}(\tau_{i,z},a_{i,z};\theta_{i})$
(22)

where $\theta_{i}$ is the network parameter.
The loss function of the VDN network is defined as

$L_{\theta}=\sum_{i}\sum_{z}(y_{i,z}^{tot}-Q_{tot}({\boldsymbol{\tau}},\mathbf{a},s_{i,z};\theta))^{2}$
(23)

$y^{tot}$ is the total target action value, which is calculated as shown in
equation (24).

$y^{tot}=r_{z}+\gamma\cdot
max_{\mathbf{a}^{\prime}}Q_{tot}({\boldsymbol{\tau}}^{\prime},\mathbf{a}^{\prime},s^{\prime}_{i,z};\theta^{\prime})$
(24)

The term $\gamma$ is the discount factor, ${\boldsymbol{\tau}}^{\prime}$,
$\mathbf{a}^{\prime}$ and $s^{\prime}_{i,z}$ are the corresponding observation
history, action and state after state $s_{i,z}$ takes action $\mathbf{a}$ and
receive reward $r_{i,z}$; $\theta^{\prime}$ is the parameter of the target
network.
Furthermore, we apply the $\epsilon-$greedy action selection to balance the
exploration-exploitation trade off, the optimal action that maximizes the
action value function is set with a probability of 1-$\epsilon$, as shown in
equation (25).

$Pr(a_{i,t,z}=a_{i,t,z}^{opt})=1-\epsilon,0<\epsilon<1$ (25)

Where $a_{i,t,z}^{opt}$ is the optimal action, and a random action is chosen
with a probability of $\epsilon$, the transition
$<s_{i,t,z},a_{i,t,z},r_{t,z},s_{i,t+1,z}>$ at time step $t$ is stored in a
replay buffer $B_{1}$. With the above information, the following price
scheduling training procedure is summarised in algorithm 1.

##  V SIMULATIONS

###  V-A simulation Setup

Figure 5: VDN accumulation reward distribution Figure 6: Q-mix accumulation
reward distribution

In our simulation set up, we adopt the Chinese time-of-use tariff for our CS 0
from [7], [4]. Furthermore, we developed the CS 1 tariff based on the original
version, with higher prices at each corresponding time period; the tariffs for
both stations are presented in Table I. The arrival time $t_{a}$ ,departure
time $t_{d}$, arrival state of charge $SOC_{a}$ and the normalized EV-CS 0
distance $\widetilde{d_{i,0}}$ (and $\widetilde{d_{i,1}}$ = 1 -
$\widetilde{d_{i,0}}$ ) all follow a normal distribution , with a
corresponding boundary limit, as shown in Table II. Moreover, the needed state
of charge $SOC^{need}$ is set to 0.8, minimum state of charge $SOC^{min}$ and
maximum state of charge $SOC^{max}$ are set to 0.2, 0.9 respectively. The
maximum charging power for NC pile is 6 kWh,with a maximum battery capacity of
$C^{NC}$ = 24 kWh, and for FC pile the maximum charging power is 30 kWh, with
a battery capacity $C^{FC}$ = 180 kWh; the random seed of 50 is selected.
For our VDN network model, the layers in the network have 64 neurons, and the
learning rate of the network is set to be 0.001. The batch size is 32, and the
buffer size is 2000.

###  V-B EV shunting evaluation

Firstly, we evaluate the effectiveness of the price scheduling strategy, i.e.
it achieves sufficient EV shunting under each case study. In this sub section,
three different scenarios are generated for the shunting evaluation, where
each scenario has varied numbers of EVs requiring charging and charging pile
sizes at each station: 1) 6 arriving EVs, CS 0 with 3 FC piles and 2 NC piles,
CS 1 with 4 FC piles and 3 NC piles; 2) 12 arriving EVs, CS 0 with 6 FC piles
and 4 NC piles, CS 1 with 8 FC piles and 6 NC piles; 3) 25 arriving EVs, CS 0
with 14 FC piles and 6 NC piles, CS 1 with 20 FC piles and 10 NC piles. The
weight coefficients for the attraction model under each scenario are
configured as

$w_{0},w_{1},w_{2}=(0.25,3.8,0.7)$ (26)

Each scenario is simulated for 1000 times, and 16 samples are intercepted for
representation, as shown in Figure 4. Under scenario 1, our price scheduling
method results in two majority arrived EV size pairs of
$(CS0,CS1)=(3,3),(4,2)$, whilst without price scheduling the common size pair
is $(5,1)$.Under scenario 2, the price scheduling results in common size pairs
of $(6,6),(7,5)$ in comparison with the baseline of $(10,2)$; and finally for
scenario 3, the result shows a difference of $(12,13)$ vs $(1,24)$. The above
results are sufficient in proving our price scheduling strategy in achieving
reasonable EV shunting.

###  V-C VDN performance evaluation

To evaluate the performance of our price scheduling strategy, we pick the VDN
model performance with the identical EV user preference model (without price
scheduling) as our baseline, and compare the corresponding accumulated rewards
under the aforementioned scenarios, depicted in Figure 5. In Table III, the
average reward at convergence and the episode of convergence are presented
under each scenario, and results of the price scheduling method, baseline are
compared. Scenario 1: Our price scheduling method converges 68 episode faster
than the baseline, while the average reward at convergence of our method is
157.537 higher than the baseline, and the accumulated rewards appear to be
more compact. Scenario 2: The price scheduling method converges 82 episodes
faster than the baseline, and the corresponding average reward is higher by
137.77; again, the accumulated rewards are virtually more compact in
comparison with the baseline results. Scenario 3: Our method converges 16
episodes faster, while achieving a significant higher average award of
-442.093, with a better overall reward compactness over the baseline.

###  V-D Q-mix performance evaluation

To further investigate the effectiveness of the proposed method, we replace
the VDN model with Q-mix network. The CS $z$ passes the agent observation
$o_{i,t,z}$ and previous action step $a_{i,t-1,z}$ to the agent network, where
the agent network is a recurrent neural network (RNN) that consists 3 layors:
multilayer perceptron (MLP), gated recurrent unit (GRU) and MLP. The agent
network outputs each value function $Q_{i,z}(\tau_{i,z},a_{i,t,z})$ into the
mixing network, where the mixing network consists of hypernetworks that take
in the state $s_{t,z}$, all the generated $Q_{i,z}(\tau_{i,z},a_{i,t,z})$
together with the weights and biases from the hypernetworks are used to
produce the total action value function
$Q_{tot}({\boldsymbol{\tau}},\mathbf{a})$.
Additionally,the local $argmax$ applied on each individual action value
function should obtain the same monotonicity effect as the global $argmax$
applied on the total action value function, in which can be expressed as

$\displaystyle argmax_{\mathbf{u}}Q_{tot}({\boldsymbol{\tau}},\mathbf{a})$
$\displaystyle=\begin{pmatrix}argmax_{u_{1,0}}Q_{1,0}(\tau_{1,0},a_{1,0})\\\
\vdots\\\ argmax_{u_{n,1}}Q_{n,1}(\tau_{n,1},a_{n,1})\end{pmatrix}$ (27)

Moreover, this expression can be further generalized to a constraint, as shown
in equation (28).

$\frac{\partial Q_{tot}}{\partial Q_{i,z}}\geq 0,\forall i\in[1..n],\forall
z\in[0..1]$ (28)

The Q-mix model adapts the same parameter settings as VDN.
The same three scenarios are implied, with the corresponding accumulated
reward performances shown in Figure 6. From Table IV, under scenario 1 the
price scheduling strategy converges 1 epoch faster over the baseline, however,
the corresponding average reward at convergence is 182.975 less than the
baseline and the accumulated rewards are shown to be less compact. Moving onto
scenario 2, the price scheduling method converges 95 episodes slower than the
baseline, but the average reward is 200.864 higher and the data scattering
appears to be more compact. For scenario 3, the price scheduling method
achieves a noticeably higher average award of -696.103 over the baseline,
similar to the simulation results for the VDN network, while achieving it with
faster episode convergence and a comparable overall reward compactness between
the two methods. The above results indicate that, under scenarios where
charging pile sizes at each station and the arrived EV sizes are sufficiently
large, the price scheduling strategy is capable to producing noticeable reward
optimization improvement over the baseline. Together with the VDN performance
results, it is sufficient to claim that under reasonably large dataset, the
proposed price scheduling strategy is capable of providing noticeable reward
optimization improvements, while achieving it with faster episode convergence.

TABLE I: Chinese Time-of-Tariff Time | CS 0 | CS 1 |
---|---|---|---
Price(CNY/kWh)
07:00 - 10:00 | 1.0044 | 1.2044 |
11:00 - 13:00 | 0.6950 | 0.7950 |
14:00 - 17:00 | 1.0044 | 1.2044 |
18:00 - 19:00 | 0.6950 | 0.7950 |
20:00 - 06:00 | 0.3946 | 0.4946 |
TABLE II: Commuting Behaviour Distributions Parameters | Distribution | Boundary |
---|---|---|---
$t_{a}$ | $\mathcal{N}(9,1^{2})$ | $7\leq t_{a}\leq 11$ |
$t_{d}$ | $\mathcal{N}(19,1^{2})$ | $17\leq t_{d}\leq 21$ |
$SOC_{a}$ | $\mathcal{N}(0.4,0.1^{2})$ | $0.2\leq SOC_{a}\leq 0.6$ |
$\widetilde{d_{i,0}}$ | $\mathcal{N}(0.5,0.3^{2})$ | $0<\widetilde{d_{i,0}}<1$ |
TABLE III: VDN Performance Results Method | Average reward | Episode of convergence |
---|---|---|---
Scenario 1
Price Schedule | -141.383 | 1170 |
Baseline | -298.92 | 1238 |
Scenario 2
Price Schedule | -63.288 | 1063 |
Baseline | -201.058 | 1145 |
Scenario 3
Price Schedule | -442.093 | 1176 |
Baseline | -547.313 | 1192 |
TABLE IV: Q-MIX Performance Results Method | Average reward | Episode of convergence |
---|---|---|---
Scenario 1
Price Schedule | -401.474 | 1214 |
Baseline | -218.499 | 1215 |
Scenario 2
Price Schedule | -112.736 | 1128 |
Baseline | -313.6 | 1223 |
Scenario 3
Price Schedule | -696.103 | 1096 |
Baseline | -1732.4 | 1142 |

##  VI CONCLUSION

To conclude, we introduced a multi-CS decentralized collaborative framework,
while employing a modified EV user preference model. Furthermore, we proposed
a price scheduling method that achieves comparable EV shunting performance.
With numerous simulations over two separate CSs and under three different
scenarios , the proposed method is proven to achieve sufficient reward
optimisation and convergence improvement results over the baseline. For
potential future research, we aim to introduce multi-neural network interface,
with the current price scheduling strategy reinforced with additional networks
to improve the capability of handling large scaled EV transportation models.
Moreover, the optimization goal can be extended, such as merging CS route
recommendation goals with price reduction, adapting further practicalities in
a real world situation.

## References

  * [1] “The uk’s plans and progress to reach net zero by 2050,” pp. 1–31, The House of Commons Library, 2023.
  * [2] “Powering up britain: Net zero growth plan,” pp. 1–126, Department for Energy Security and Net zero, 2023.
  * [3] H. Song, C.-C. Liu, J. Lawarrée, and R. W. Dahlgren, “Optimal electricity supply bidding by markov decision process,” IEEE transactions on power systems, vol. 15, no. 2, pp. 618–624, 2000.
  * [4] H. Li, G. Li, T. T. Lie, X. Li, K. Wang, B. Han, and J. Xu, “Constrained large-scale real-time ev scheduling based on recurrent deep reinforcement learning,” International Journal of Electrical Power & Energy Systems, vol. 144, p. 108603, 2023.
  * [5] Z. Wen, D. O’Neill, and H. Maei, “Optimal demand response using device-based reinforcement learning,” IEEE Transactions on Smart Grid, vol. 6, no. 5, pp. 2312–2324, 2015.
  * [6] S. Iqbal and F. Sha, “Actor-attention-critic for multi-agent reinforcement learning,” in International conference on machine learning, pp. 2961–2970, PMLR, 2019.
  * [7] H. Li, B. Han, G. Li, K. Wang, J. Xu, and M. W. Khan, “Decentralized collaborative optimal scheduling for ev charging stations based on multi-agent reinforcement learning,” IET Generation, Transmission & Distribution, 2023.
  * [8] L. Zou, M. S. Munir, Y. K. Tun, S. Kang, and C. S. Hong, “Intelligent ev charging for urban prosumer communities: An auction and multi-agent deep reinforcement learning approach,” IEEE Transactions on Network and Service Management, vol. 19, no. 4, pp. 4384–4407, 2022.
  * [9] T. Sousa, H. Morais, Z. Vale, and R. Castro, “A multi-objective optimization of the active and reactive resource scheduling at a distribution level in a smart grid context,” Energy, vol. 85, pp. 236–250, 2015.
  * [10] T. Qian, C. Shao, X. Wang, and M. Shahidehpour, “Deep reinforcement learning for ev charging navigation by coordinating smart grid and intelligent transportation system,” IEEE transactions on smart grid, vol. 11, no. 2, pp. 1714–1723, 2019.
  * [11] Q. Xing, Y. Xu, Z. Chen, Z. Zhang, and Z. Shi, “A graph reinforcement learning-based decision-making platform for real-time charging navigation of urban electric vehicles,” IEEE Transactions on Industrial Informatics, vol. 19, no. 3, pp. 3284–3295, 2022.
  * [12] B. Wang, B. Hu, C. Qiu, P. Chu, and R. Gadh, “Ev charging algorithm implementation with user price preference,” in 2015 IEEE Power & Energy Society Innovative Smart Grid Technologies Conference (ISGT), pp. 1–5, IEEE, 2015.
  * [13] J. Liu, S. Wang, and X. Tang, “Pricing and charging scheduling for cooperative electric vehicle charging stations via deep reinforcement learning,” in 2022 IEEE International Conference on Communications, Control, and Computing Technologies for Smart Grids (SmartGridComm), pp. 212–217, IEEE, 2022.
  * [14] P. Sunehag, G. Lever, A. Gruslys, W. M. Czarnecki, V. Zambaldi, M. Jaderberg, M. Lanctot, N. Sonnerat, J. Z. Leibo, K. Tuyls, et al., “Value-decomposition networks for cooperative multi-agent learning,” arXiv preprint arXiv:1706.05296, 2017.
  * [15] T. Rashid, M. Samvelyan, C. S. De Witt, G. Farquhar, J. Foerster, and S. Whiteson, “Monotonic value function factorisation for deep multi-agent reinforcement learning,” The Journal of Machine Learning Research, vol. 21, no. 1, pp. 7234–7284, 2020.
  * [16] Y. Zhao, Y. Guo, Q. Guo, H. Zhang, and H. Sun, “Deployment of the electric vehicle charging station considering existing competitors,” IEEE Transactions on smart grid, vol. 11, no. 5, pp. 4236–4248, 2020.