# Superconductivity in Nb: from the impact of temperature, Cooper-pairing to
dimensionality

Uriel A. Aceves Rodriguez1,2 Filipe Guimarães3 Samir Lounis1,2 s.lounis@fz-
juelich.de 1Peter Grünberg Institut and Institute for Advanced Simulation,
Forschungszentrum Jülich & JARA, 52425 Jülich, Germany
2Faculty of Physics & CENIDE, University of Duisburg-Essen, 47053 Duisburg,
Germany
3Jülich Supercomputing Centre, Forschungszentrum Jülich & JARA, 52425 Jülich,
Germany

###### Abstract

The ability to simulate realistically the electronic structure of
superconducting materials is important to understand and predict various
properties emerging in both the superconducting topological and spintronics
realms. We introduce a tight-binding implementation of the Bogoliubov-de
Gennes method, parameterized from density functional theory, which we utilize
to explore the bulk and thin films of Nb, known to host a significant
superconducting gap. The latter is useful for various applications such as the
exploration of trivial and topological in-gap states. Here, we focus on the
simulations aspects of superconductivity and study the impact of temperature,
Cooper-pair coupling and dimensionality on the size of the superconducting
order parameter.

##  I Introduction

Superconductivity is a basic fundamental phenomenon in solid-state physics.
Although being discovered more than a century ago, it still generates a
plethora of research activities. Recently, its importance flourished again in
the context of topological superconductivity[1, 2, 3], which is an essential
field of research for the realization of quantum computing. Certainly, the
intricate interplay of standard spintronics and quantum information concepts
can give rise to unforeseen breakthroughs in basic research and in future
technologies for storing, computing and transmitting information. This calls
for the development of theoretical frameworks based on a realistic description
of the electronic structure of superconducting materials interfaced with
magnetic systems in order to understand and predict various related phenomena.

Interesting topics emerged in superconducting spintronics [4, 5, 6, 7], where
the interplay of Cooper pairs and other electronic degrees of freedom (such as
spin, charge and orbital) give rise to fascinating phenomena, ranging from
supercurrents, triplet superconductivity, and torques. Simultaneously, the
potential design of topological superconductivity and the creation of Majorana
states[8, 9, 10] is driving a lot of research in the context of topological
qubits [11].

The microscopic theory to describe conventional superconductivity goes back to
the work of John Bardeen, Leon Cooper, and Robert Schrieffer (BCS)[12]. In
this context, the Bogoliubov–de Gennes (BdG) method[13, 14, 15, 16] is an
elegant mean-field approximation that relies upon Bogoliubov–Valatin
transformations that take the Hamiltonian from a particle space into a
particle-hole one—a framework that has been frequently used, see for example
Refs. 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27. Nowadays, several methods
based on a realistic description of the electronic structure of materials are
capable of treating fundamental aspects related to superconductivity. For
instance, powerful theoretical frameworks have been proposed and utilized
within density functional theory (DFT) combined with the BdG formalism[28, 29,
30, 31, 32, 33], or by proposing an extension of DFT to account for
superconductivity [34, 35, 36, 37]. Tight-binding schemes are very appealing
since they allow a versatile description of superconductivity-related physics
with a potential control of various levels of approximations, which could
enable the treatment of problems not attainable with conventional DFT.

Among the superconducting materials, Nb stands out and is often utilized as a
key component to investigate superconductivity-related physics [38, 39]. The
surface of Nb(110) became recently the standard substrate [40], where Yu-
Shiba-Rusinov in-gap states[9, 8, 10] and the topological Majorana states[8,
9, 10] are explored in adatoms[41, 42, 29, 43], at the edges of magnetic
wires[44, 24, 26, 27], or films [45, 46, 47, 48].

The goal of this article is to introduce our recently developed method that
allows the treatment of superconductivity on the basis of a tight-binding
approach as implemented in TITAN [49, 50, 51, 52]. The framework enables a
realistic description of the electronic structure, with parameters obtained
from DFT [53]. Examples of applications of TITAN ranges from dynamical
magnetic responses[54], dynamical transport and torques[50, 52], magnetic
damping [51], as well as ultrafast magnetization dynamics triggered by laser
pulses [55, 56]. In the context of superconductivity, the method was already
introduced in Ref. 57, with a focus on magnetic exchange interactions at
superconducting-magnetic interfaces. Here, we pay attention to Nb in bulk and
thin films with the aim of highlighting fundamental aspects in the simulations
procedure. Aspects related to the choice of the electron-phonon coupling to
open the appropriate band gap are discussed together with the impact of
temperature and self-consistency on the underlying physics. Dimensionality
affects superconductivity, as we unveil, which can be utilized to tune the
superconducting gap for different applications in modern investigations of
superconductivity-related problems.

The article is organized as follows. We first introduce the method in Sec. II,
which is based on multi-orbital tight-binding theory. Then we present in Sec.
III the results of our simulations in Nb bulk, (001) and (110) monolayers,
which are followed by Nb(110) thin films of different thicknesses. Finally, we
present our concluding remarks in Sec. IV.

##  II Tight-binding and the Bogoliubov–de Gennes method

To simulate the superconducting properties of Nb across dimensions, we
introduce here the tight-binding methodology together with the BdG method that
we utilized, as implemented in our code TITAN [49, 51, 52, 55, 56, 57]. The
hopping parameters were obtained from first-principles calculations [53]. More
details on superconducting related aspects are elaborated in Ref. 57.

The fundamental tight-binding Hamiltonian reads:

$\displaystyle
H_{S}=\frac{1}{N}\left\\{\sum_{\alpha\beta,\sigma\eta,\mu\nu,\bm{k}}H_{\alpha\beta,\sigma\eta}^{\mu\nu}(\bm{k})c^{\dagger}_{\alpha\mu\sigma}(\bm{k})c_{\beta\nu\eta}(\bm{k})-\sum_{\alpha,\mu\nu,\bm{k}\bm{k}^{\prime}}\lambda_{\alpha\mu}c^{\dagger}_{\alpha\mu\uparrow}(\bm{k})c^{\dagger}_{\alpha\mu\downarrow}(-\bm{k})c_{\alpha\mu\downarrow}(-\bm{k}^{\prime})c_{\alpha\mu\uparrow}(\bm{k}^{\prime})\right\\}\,,$
(1)

where $c^{\dagger}_{\alpha\mu\sigma}(\bm{k})$ and $c_{\beta\nu\eta}(\bm{k})$
consist of the creation and annihilation operators associated to electrons
having a wave vector $\bm{k}$ and spin $\sigma$ in the orbitals $\mu$ and spin
$\eta$ in the orbitals $\nu$ , respectively. $\alpha$ and $\beta$ are specific
atoms in the bulk unit cell, or layers in the film geometry. $\bm{k}$ is a
reciprocal vector, while $N$ is the number of wave vectors in the three-
dimensional (for bulk) or two-dimensional (for films) Brillouin zone. To
enable the formation of Cooper pairs, which ultimately gives rise to
superconductivity, we introduce $\lambda_{\alpha\mu}$ in the second term of
the previous equation, as induced by the electron-phonon coupling.

In the mean-field approximation, Eq. (1) simplifies to

$\displaystyle H_{S}^{\text{MF}}=$
$\displaystyle\frac{1}{N}\sum_{\bm{k}}\left\\{\sum_{\alpha\beta,\sigma\eta,\mu\nu}H_{\alpha\beta,\sigma\eta}^{\mu\nu}(\bm{k})c^{\dagger}_{\alpha\mu\sigma}(\bm{k})c_{\beta\nu\eta}(\bm{k})\right.$
(2)
$\displaystyle-\left.\sum_{\alpha,\mu}\left(\Delta_{\alpha\mu}^{*}c_{\alpha\mu\downarrow}(-\bm{k})c_{\alpha\mu\uparrow}(\bm{k})+\Delta_{\alpha\mu}c_{\alpha\mu\uparrow}^{\dagger}(\bm{k})c_{\alpha\mu\downarrow}^{\dagger}(-\bm{k})\right)\right\\},$

with

$\Delta_{\alpha\mu}=\lambda_{\alpha\mu}\frac{1}{N}\sum_{\bm{k}}\langle
c_{\alpha\mu\downarrow}(-\bm{k})c_{\alpha\mu\uparrow}(\bm{k})\rangle,\quad\Delta_{\alpha\mu}^{*}=\lambda_{\alpha\mu}\frac{1}{N}\sum_{\bm{k}}\langle
c_{\alpha\mu\uparrow}^{\dagger}(\bm{k})c_{\alpha\mu\downarrow}^{\dagger}(-\bm{k})\rangle,$
(3)

$\Delta_{\alpha\mu}$ is known as the superconducting gap parameter, which
corresponds to the energy required to scatter Cooper pairs [58].

The term $H_{\alpha\beta,\sigma\eta}^{\mu\nu}(\bm{k})$ corresponds to the
regular Hamiltonian where superconductivity is not accounted for:

$H^{\mu\nu}_{\alpha\beta,\sigma\eta}(\bm{k})=H_{\alpha\beta}^{0\mu\nu}(\bm{k})\sigma^{0}+\bm{\sigma}\cdot\hat{\bm{e}}_{\alpha}B_{\alpha}^{[\text{xc}]\mu\nu}(\bm{k})\delta_{\alpha\beta}+\bm{\sigma}\cdot\bm{B}_{\alpha}^{[\text{soc}]\mu\nu}(\bm{k})\delta_{\alpha\beta}\,,$
(4)

with $H_{\alpha\beta}^{0\mu\nu}$ being the spin-independent tight-binding
term, the second term is responsible for the intra-atomic exchange interaction
(originating from a Hubbard-like contribution[50, 52]), which can lead to
magnetism, and the third term corresponds to the spin-orbit interaction.

The BdG methodology consists in utilizing the Bogoliubov–Valatin
transformation[16] to rewrite Eq. (2) in the electron-hole representation
instead of the electron one, which leads to the BdG Hamiltonian:

$H_{\text{BdG}}^{\alpha\beta,\mu\nu}(\bm{k})=\begin{pmatrix}H_{\alpha\beta,\uparrow\uparrow}^{\mu\nu}(\bm{k})-E_{F}&H_{\alpha\beta,\uparrow\downarrow}^{\mu\nu}(\bm{k})&0&-\Delta_{\alpha\mu}\mathbb{I}\\\
H_{\alpha\beta,\downarrow\uparrow}^{\mu\nu}(\bm{k})&H_{\alpha\beta,\downarrow\downarrow}^{\mu\nu}(\bm{k})-E_{F}&\Delta_{\alpha\mu}\mathbb{I}&0\\\
0&\Delta_{\alpha\mu}^{*}\mathbb{I}&-H_{\alpha\beta,\uparrow\uparrow}^{\mu\nu*}(-\bm{k})+E_{F}&-H_{\alpha\beta,\uparrow\downarrow}^{\mu\nu*}(-\bm{k})\\\
-\Delta_{\alpha\mu}^{*}\mathbb{I}&0&-H_{\alpha\beta,\downarrow\uparrow}^{\mu\nu*}(-\bm{k})&-H_{\alpha\beta,\downarrow\downarrow}^{\mu\nu*}(-\bm{k})+E_{F}\end{pmatrix},$
(5)

with $E_{F}$ being the Fermi energy while the associated eigenvalue problem to
solve reads:

$\sum_{\beta\mu}H_{\text{BdG}}^{\alpha\beta,\mu\nu}(\bm{k})\phi_{\beta\mu}(\bm{k})=E_{n}(\bm{k})\phi_{\alpha\nu}(\bm{k}).$
(6)

The transformation is canonical and enabled by rewriting

$c_{\alpha\mu\sigma}(\bm{k})=\sum_{n}^{{}^{\prime}}u_{\alpha\sigma}^{n}(\bm{k})\gamma_{n}+v_{\alpha\sigma}^{n*}(\bm{k})\gamma_{n}^{\dagger},\quad
c_{\alpha\mu\sigma}^{\dagger}(\bm{k})=\sum_{n}^{{}^{\prime}}u_{\alpha\sigma}^{n*}(\bm{k})\gamma_{n}^{\dagger}+v_{\alpha\sigma}^{n}(\bm{k})\gamma_{n}^{\dagger},$
(7)

with

$\\{\gamma_{n},\gamma_{m}\\}=\\{\gamma_{n}^{\dagger},\gamma_{m}^{\dagger}\\}=0,\quad\\{\gamma_{n}^{\dagger},\gamma_{m}\\}=\delta_{nm},$
(8)

where the prime indicates that the sums run only over the states with positive
energy[16, 59]. This restriction in the sum is done to counteract the doubling
of the degrees of freedom originated from the change of basis.

The eigenvector of Eq. (6) is

$\phi_{i\nu}(\bm{k})=\begin{pmatrix}u_{\alpha\nu\uparrow}(\bm{k})\\\
u_{\alpha\nu\downarrow}(\bm{k})\\\ v_{\alpha\nu\uparrow}(\bm{k})\\\
v_{\alpha\nu\downarrow}(\bm{k})\end{pmatrix}.$ (9)

Our numerical procedure to evaluate the various properties of a given material
is based on the self-consistency of the charge density, which is also pursued
once the electron-phonon coupling is included. In practice, we start the
simulations considering the normal metallic phase then we incorporate the
Cooper pairing and repeat the computational runs until convergence.

##  III Results

In this section we explore the results of simulations done for various Nb
structures. Our choice on Nb is motivated by its relatively large
superconducting gap, which facilitates the exploration of in-gap states.
Simulating low temperatures phenomena requires finer computational meshes,
thus for the first investigations we choose high temperatures and large values
for the superconducting coupling $\lambda$. We start with bulk Nb, where we
recover a superconducting gap opening in the electronic structure around the
Fermi energy. Then we explore the case of thin films down to the single
monolayer along the [001] or [110] directions. We monitor how the
superconducting gap is affected by the self-consistency of the electronic
structure calculations.

###  III.1 Bulk superconducting Nb

Nb has a superconducting critical temperature of
$T_{C}=$9.3\text{\,}\mathrm{K}$$[60], which leads to a significant
superconducting gap $2\Delta=$3.8\text{\,}\mathrm{meV}$$. A rough
approximation relates $T_{c}$ with the size of the gap—at T=0 K—through [58]

$\Delta(T=0)=1.764\,k_{b}T_{c}\,.$ (10)

Large values of $\Delta$ are desirable when investigating states that appear
inside the gap. In our simulation procedure, we start with bulk Nb at high
temperatures. In TITAN, the thermal effects enter through the parameter
$\eta=k_{b}T/\pi$, which broadens the Fermi-Dirac function; in this case, we
will use $k_{b}T=$43\text{\,}\mathrm{meV}$$ ($\eta=$1\text{\,}\mathrm{mRy}$$).
For a reference, this is higher than room temperature, which is at
$26\text{\,}\mathrm{meV}$, but nevertheless, considerably lower than the
melting point of Nb $237\text{\,}\mathrm{meV}$. Nb has a body-centered cubic
crystal structure. We make use of the Slater-Koster (SK) parameters from Ref.
53, the spin-orbit coupling strength comes from Ref. 61. In Fig. 1(a), we show
the density of states of both spin channels (left and right panels), and in
the middle the band structure of bulk Nb. The charge occupancy per orbital is
$\langle n_{s}\rangle=0.70\,,\langle n_{p}\rangle=0.68\,,\langle
n_{d}\rangle=3.61.$ The density of states (DOS) and the band structure agree
with the results displayed in Ref. 53.

Figure 1: Electronic structure of bulk Nb with and without superconductivity.
Three cases are illustrated: a without the BdG formalism only the electron
part is shown; b with the BdG formalism but setting the Cooper pair coupling
$\lambda$ to zero; c with BdG and $\lambda=$5.44\text{\,}\mathrm{eV}$$ in all
orbitals. The red colored line evinces the location of the gap as predicted by
the self-consistent value $\Delta=$0.497\text{\,}\mathrm{eV}$$ and thermal
broadening $\eta=$43\text{\,}\mathrm{meV}$$. The band structure is shown in
the middle. To the right (left), the DOS of the majority (minority) spin
channel is presented. The legend box depicts the colors used to signal the
total, the band resolved DOS and their electron versus hole parts, i.e. the
contributions from respectively the $u$ and $v$ components of the
eigenvectors.

Repeating the calculations by solving the BdG equations but setting
$\lambda=0$ on all orbitals, we naturally recover the metallic band-structure
shown in Fig. 1(b), where both the electron and hole parts are illustrated.
Both parts are symmetric with respect to the Fermi level. This holds true for
the DOS as well. We expect that by allowing for finite coupling $\lambda$,
degeneracies occurring at crossing points will be lifted.

In subsection III.2, we describe a rather optimized search method to tune this
parameter. From our explorations, we found that: (i) large values of $\lambda$
do not always mean large gaps and (ii) the size of the gap is not a linear
function of $\lambda$. From tabulated values of the electron-phonon
$\lambda_{\text{e-p}}$ coupling constant we can evaluate the coupling
parameter $\lambda$ utilizing the following relation [62]

$\lambda=\lambda_{\text{e-p}}D(E_{F})\,,$ (11)

where $D(E_{F})$ is the density of states at the Fermi energy. For example,
for bulk Nb with the data from Refs. 63 we obtain
$\lambda\approx$1.3\text{\,}\mathrm{eV}$$. In Ref. 62, the electron-phonon
coupling constant is calculated for slabs of Nb(100) and Nb(110) of 3, 6, 9,
12, and 15 layers. For the Nb(110), the results for the system of 9 layers—the
slab size closer to the one we will study later—lead to coupling parameters of
$\lambda\approx$5.28\text{\,}\mathrm{eV}$$ at the surfaces and
$\lambda\approx$1.3\text{\,}\mathrm{eV}$$ in the middle layer. When plotting
the curve associated to $\Delta$ vs. $\lambda$, we obtain a quadratic
function, in agreement with Ref. 64, where the self-consistent
Korringa–Kohn–Rostoker method is used. There, the authors found that in order
to obtain $\Delta\approx$1\text{\,}\mathrm{meV}$$ the coupling parameter must
be $\lambda\approx$1.1\text{\,}\mathrm{eV}$$ for bulk Nb.

As a first example of our investigations, we set
$\lambda_{\mu}=$5.44\text{\,}\mathrm{eV}$$. The occupancies are $\langle
n_{s}\rangle=0.70317397\,,\langle n_{p}\rangle=0.66860504\,,\langle
n_{d}\rangle=3.62822099$, which experience a minor change of the order of
$10^{-2}$ with respect to the case where $\lambda$ is set to zero. The self-
consistent value of the superconducting gap parameter is
$\Delta=$0.497\text{\,}\mathrm{eV}$$, which is obviously gigantic. The
resulting band structure and LDOS plots are illustrated in Fig. 1(c). The red
dotted lines delimit the range $E\in[-\Delta,\Delta]$. As a sanity check, it
is reassuring that the gap in the band structure matches the self-consistent
output value of $\Delta$. The density of states does not show a perfect gap,
but decays at the gap region. This effect is due to the thermal broadening
induced artificially, which is in this case set to $43\text{\,}\mathrm{meV}$.

In subsections III.2 and III.4, we will systematically explore the impact of
temperature and coupling strengths, searching for an optimal $\lambda$ that
yields the experimental value of $\Delta$.

###  III.2 Superconducting gap tuning: single bcc (001) Nb monolayer

Explorations to find the appropriate coupling $\lambda$ that results in the
desired gap parameter $\Delta$ involve systematic repeated simulations with
different parameters. In this section, we show the typical curves for $\Delta$
as a function of the temperature (via the broadening parameter) and the Cooper
pair coupling parameter. To exemplify, we use a 2D system, namely a (001)
monolayer of bcc Nb. As in our previous examples, we set the same $\lambda$
for all orbitals, with eight values between $4\text{\,}$ and
$9\text{\,}\mathrm{eV}$. We consider forty different temperature points
between $300\text{\,}\mathrm{K}$ and $14\,000\text{\,}\mathrm{K}$. Simulations
at lower temperatures are computationally highly expensive since they require
an increased number of k-points to be able to capture the narrow peaky
structures in the Brillouin zone with lower broadening. All the simulations in
this section use the same number of k-points, $\approx 2\times 10^{6}$.

Figure 2: Dependence of the superconducting gap of Nb(001) monolayer on
various parameters. a Heat map of the superconducting gap as function of the
coupling $\lambda$ and the thermal broadening parameter. All curves drop to
zero after reaching their corresponding critical temperature. b Another
version of a in order to better follow how the superconducting gap evolves
against temperature. The bands in color reflect an rough extrapolation to the
critical temperature. c Evolution of the superconducting gap across the self-
consistency for several coupling strengths $\lambda$. d Relative change of the
superconducting gap when comparing the values obtained after one single
iteration to those after self-consistency
($(\Delta_{\text{initial}}-\Delta_{\text{final}})/\Delta_{\text{initial}}$)
for different values of $\lambda$.

Fig. 2(a) displays a heatmap with the results of the gap parameter for the
different temperatures and couplings. We can see that the gap parameter drops
to zero after reaching the corresponding critical temperature for every line.
Directing our attention for fixed temperatures (i.e. the “columns” in the
figure), it is clear that there is also a critical $\lambda_{c}$ below which
$\Delta$ is zero. Throughout the heatmap, $\Delta$ ranges from 0.0 to
$1.76\text{\,}\mathrm{eV}$. We can analyze the same data from a different
perspective, as shown in Fig. 2(b). There we see that the $\Delta$ vs. $T$
curve has indeed the shape predicted by the BCS theory. For sufficiently low
temperatures, $\Delta$ approaches a plateau. Around the critical temperature
$T_{c}$, the self-consistency procedure may become unstable and difficult to
achieve a precise value, which lead some curves to miss points in that range.
We can proceed to an extrapolation to evaluate the missing points on the basis
of the asymptotic behavior of $\Delta$ around the critical temperature [58]

$\Delta(T)\approx 1.74\Delta(0)\sqrt{1-\frac{T}{T_{c}}}\,.$ (12)

We used the precedent equation to go from the data in Fig. 2(b) to the one in
Fig. 2(a). The BCS theory predicts that, at zero Kelvin, the relation between
$\Delta$ and $T_{c}$ is the one in Eq. (10). Experimentally the reported
values range from $2\Delta=3.2k_{b}T_{c}$ to $2\Delta=4.6k_{b}T_{c}$, we
signaled this range with the colored bands in Fig. 2(b). This plot establishes
the aforementioned statement regarding the non-linearity of $\Delta$ with
respect to $\lambda$.

For a wide range of simulations, it is desirable and recommendable to place
the superconducting system well below $T_{c}$, at a point in the plateau of
the $\Delta$ vs. $T$ curve. The superconducting state is frailer at the
vicinity of the critical temperature, and internal changes in the material can
perturb the properties of interest.

###  III.3 Renormalization after self-consistency: single bcc (001) Nb
monolayer

A common problem in simulating superconductors interfaced with normal
conductors is the artificially increased size of the gap in order to
facilitate the exploration of emerging in-gap states and it is rather common
to proceed to single-iteration simulations, which assumes the electronic
structure is not significantly altered in comparison to the case, where
instead one proceeds to self-consistency. The importance of self-consistency
was highlighted in Ref. 65. In this subsection, we address this aspect of
simulation of superconducting materials by comparing the size of the gap
obtained after one iteration (starting from the converged electronic structure
of the metallic phase) with respect to the one obtained after full self-
consistency. In Fig. 2(c), we single out the value of the superconducting gap
parameter through several simulations with different coupling parameters at a
temperature of $43\text{\,}\mathrm{meV}$. Our algorithm, based on the Powell
method [66] can halt the simulations when no improvement is detected in the
last five steps. All shown cases converged without a problem.

In Fig. 2(d) we portray the percentage change of the gap parameter
($(\Delta_{\text{initial}}-\Delta_{\text{final}})/\Delta_{\text{initial}}$)
between the first- and last-iteration calculated $\Delta$ for different
$\lambda$. We note that the change is bigger for smaller values of $\lambda$.

So far, we have been discussing the numerics of the simulations assuming both
large superconducting gaps and large electronic temperatures. In section
III.4, we drop $T$ to realistic values associated to a superconducting Nb
layer. Additionally, we explore the impact on the orbital-dependent gaps
$\Delta_{\mu}$.

###  III.4 Superconducting Nb (110) monolayer

In this subsection we consider the cryogenic regime, which calls for
computationally heavier simulations. As we hinted in previous subsections,
sinking the temperature forces us to increase the number of k-points for each
simulation to be able to capture the narrow states. Since these calculations
will be more intensive, we will leave aside the testing structures we were
using previously—bulk Nb and Nb (100) monolayers—and move to our main
interest, which is Nb (110) films, heavily explored in various recent
experiments.

We assume a thermal broadening of $4.27\text{\,}\mathrm{meV}$ (equivalent to a
temperature of approximately $4.56\text{\,}\mathrm{K}$), shallow enough to be
in the plateau of the $\Delta$ vs. $T$ curve for all the $\lambda$ we used. We
consider a non-uniform grid with 57 points for $\lambda$, ranging from
$0.68\text{\,}\mathrm{eV}$ to $13.6\text{\,}\mathrm{eV}$, with a higher
resolution for lower energies. In Fig. 3(a), we display the results of the
orbital-dependent $\Delta_{\mu}$ for the entire range. Obviously, some
orbitals are more responsive to the coupling $\lambda$. This is affected by
the strength of the latter coupling. We notice an increment in $\Delta$ as we
enlarge the coupling strength. $\lambda>$2.5\text{\,}\mathrm{eV}$$ induces
gaps of the order of $1\text{\times}{10}^{2}\text{\,}\mathrm{eV}$ or larger,
high enough to be considered unrealistic.

Figure 3: Orbital-resolved superconducting gap for Nb(110) monolayer. a
$\Delta_{\mu}$ grows considerably for all the orbital except $p_{z}$. b Two
close ups to the region where the material transitions from the non-
superconducting to the superconducting phase. Simulations performed at
$4.56\text{\,}\mathrm{K}$.

In Fig. 3(b), we focus only on the low energy range
$\lambda\in[$0.5\text{\,}\mathrm{eV}$,$0.75\text{\,}\mathrm{eV}$]$, where the
growth of $\Delta(\lambda)$ is more controlled. In the inset, we also restrict
the values on the $y$-axis from $0$ to $0.8\text{\,}\mathrm{eV}$ to resolve
better the smaller gap parameter for some of the orbitals. We observe that the
$p_{z}$ component is still quite small; this is a consequence of the shape of
this orbital and the geometry of the monolayer. The orientation of the $p_{z}$
orbital is out of the plane, and this means its in-plane component is zero,
and it is hard for superconducting electrons in this orbital to interact with
others.

The gap is opened around the Fermi surface when the $u$ and $v$ states mix due
to the off-diagonal components of the Hamiltonian; this means that the gap
must come from the bands that cross $E_{F}$. Hence, the gap of the system
comes from the minimum value of $\Delta_{\mu}$ among orbitals whose bands
traverse the Fermi level. Once we profile these bands, we can execute any
searching algorithm—e.g., binary search, or nested intervals methods—until we
find an appropriate $\lambda$. As a rule of thumb, each significant figure in
$\Delta$ requires two on $\lambda$; In other words, having $\Delta$ up to
three significant digits typically requires knowing $\lambda$ to six. Refining
can be very expensive. In our case, $\lambda=$1.37\text{\,}\mathrm{eV}$$
induces $\Delta=$1.4\text{\,}\mathrm{meV}$$.

In the next subsection, we explore the case of thin films of Nb (110) at
$4.56\text{\,}\mathrm{K}$ and analyse the evolution of the superconducting gap
parameter across the material.

###  III.5 Superconducting slab

It is difficult to predict if the coupling strength $\lambda$ that works well
for one given system still opens a reasonable gap for a different one. We
started with the value found in the previous section,
$\lambda=$1.37\text{\,}\mathrm{eV}$$, which opened a gap
$\Delta=$1.4\text{\,}\mathrm{meV}$$ for a single monolayer. It turns out,
however, that as soon as the Nb film thickness is larger than 5 layers, a
large coupling is required to open a gap. In the remained study, we assume
$\lambda=$2.72\text{\,}\mathrm{eV}$$ at a temperature of
$4.56\text{\,}\mathrm{K}$.

In Fig. 4(a), we show the magnitude of the gap in the central layer of the
films as a function of the slab thickness. The orange line corresponds to the
gap parameter of the bulk Nb for $\lambda=$2.72\text{\,}\mathrm{eV}$$ at the
same temperature. One notices that the thicker the film gets, the closer the
gap is to the one obtained for the bulk, a behavior comparable to what has
been reported in Ref. 67.

Figure 4: Superconducting gap for Nb(110) films. a Superconducting gap at the
centers of Nb(110) films of different thicknesses. b Layer-resolved
superconducting gap for a nine-layer-thick Nb(110) film. The orange line
represents the gap ($46.9\text{\,}\mathrm{meV}$) associated to bulk Nb
assuming $\lambda=$2.72\text{\,}\mathrm{eV}$$ and a temperature of
$4.56\text{\,}\mathrm{K}$.

The superconducting gap is not homogeneous across the Nb films. An example is
shown in Fig. 4(b) for a film with a thickness of 9 layers. The gap in the
middle of the film is larger than at the surface and it seems to oscillate as
a function of position, as we would expect for wave-like particles in a
confined system. A lower gap at the surface is expected since the electrons
there have fewer opportunities to couple in Cooper-pairs due to their lack of
neighbors from one side. We must also consider surface effects; for example,
we have around 1.5 electrons less than the inner parts. Electrons move from
the surfaces to the middle layers. Interestingly some of the layers in the
film can have a larger gap than the bulk one. This indicates that the
dimensionality of a material, can be utilized to control the superconducting
gap due to confinement effects.

##  IV Conclusions

After briefly introducing the tight-binding method that we developed to tackle
superconductivity within the Bogoliubov-de Gennes formalism, we presented the
results of our simulations on a material that is currently heavily utilized to
explore the physics of superconductivity. We addressed the case of bulk and
thin films of Nb along different directions. A major issue in this field from
the computational point of view is to choose the right Cooper pair coupling
that opens the desired superconducting gap since the calculations are
extremely expensive to resolve gaps of a few $\text{\,}\mathrm{meV}$ or even
smaller. We analysed how the superconducting order parameter correlates with
both the Cooper pair coupling and temperature, and explored the impact of
self-consistency on the emerging orbital-dependent superconducting gap. In the
thin film geometry, we unveiled the possibility of engineering the magnitude
of the gap via confinement effects. Indeed, the superconducting gap is found
to be layer-dependent, either smaller or larger than the values found in the
bulk phase.

###### Acknowledgements.

The authors acknowledge funding provided by the Priority Programmes SPP 2137
“Skyrmionics” (Projects LO 1659/8-1) of the Deutsche Forschungsgemeinschaft
(DFG). We acknowledge the computing time provided through JARA on the
supercomputers JURECA[68]. Simulations were also performed with computing
resources granted by RWTH Aachen University under project jara0189 and
p0020362.

## References

  * Lahtinen and Pachos [2017] V. Lahtinen and J. K. Pachos, A Short Introduction to Topological Quantum Computation, SciPost Phys. 3, 021 (2017).
  * Mishra  _et al._ [2021] A. Mishra, P. Simon, T. Hyart, and M. Trif, Yu-Shiba-Rusinov qubit, PRX Quantum 2, 040347 (2021).
  * Pavešić and Žitko [2022] L. Pavešić and R. Žitko, Qubit based on spin-singlet Yu-Shiba-Rusinov states, Phys. Rev. B 105, 075129 (2022).
  * Linder and Robinson [2015] J. Linder and J. W. A. Robinson, Superconducting spintronics, Nature Physics 11, 307 (2015).
  * Eschrig [2015] M. Eschrig, Spin-polarized supercurrents for spintronics: a review of current progress, Reports on Progress in Physics 78, 104501 (2015).
  * Buzdin [2005] A. I. Buzdin, Proximity effects in superconductor-ferromagnet heterostructures, Rev. Mod. Phys. 77, 935 (2005).
  * Bergeret  _et al._ [2018] F. S. Bergeret, M. Silaev, P. Virtanen, and T. T. Heikkilä, Colloquium: Nonequilibrium effects in superconductors with a spin-splitting field, Rev. Mod. Phys. 90, 041001 (2018).
  * Beenakker [2013] C. Beenakker, Search for Majorana fermions in superconductors, Annual Review of Condensed Matter Physics 4, 113 (2013), https://doi.org/10.1146/annurev-conmatphys-030212-184337 .
  * Kitaev [2001] A. Y. Kitaev, Unpaired Majorana fermions in quantum wires, Physics-Uspekhi 44, 131 (2001).
  * Tanaka  _et al._ [2009] Y. Tanaka, T. Yokoyama, and N. Nagaosa, Manipulation of the majorana fermion, andreev reflection, and josephson current on topological insulators, Phys. Rev. Lett. 103, 107002 (2009).
  * Aguado and Kouwenhoven [2020] R. Aguado and L. P. Kouwenhoven, Majorana qubits for topological quantum computing, Physics Today 73, 44 (2020), https://pubs.aip.org/physicstoday/article-pdf/73/6/44/10124104/44_1_online.pdf .
  * Bardeen  _et al._ [1957] J. Bardeen, L. N. Cooper, and J. R. Schrieffer, Microscopic theory of superconductivity, Phys. Rev. 106, 162 (1957).
  * Bogoliubov [1958] N. Bogoliubov, A New Method in the Theory of Superconductivity. I, Journal of Experimental and Theoretical Physics 7, 41 (1958).
  * De Gennes and Pincus [1999a] P. De Gennes and P. Pincus,  _Superconductivity Of Metals And Alloys_, Advanced Books Classics (Avalon Publishing, 1999).
  * Alaberdin  _et al._ [1996] E. R. Alaberdin, A. A. Vichorev, A. M. Savchenko, and B. I. Sadovnikov, On Bogoliubov's method in superconductivity, Theoretical and Mathematical Physics 107, 523 (1996).
  * Zhu [2016] J. Zhu,  _Bogoliubov-de Gennes Method and Its Applications_ , Lecture Notes in Physics (Springer International Publishing, 2016).
  * Han [2009] Q. Han, A method of studying the bogoliubov–de gennes equations for the superconducting vortex lattice state, Journal of Physics: Condensed Matter 22, 035702 (2009).
  * Saunderson  _et al._ [2020a] T. G. Saunderson, Z. Győrgypál, J. F. Annett, G. Csire, B. Újfalussy, and M. Gradhand, Real-space multiple scattering theory for superconductors with impurities, Phys. Rev. B 102, 245106 (2020a).
  * Pellegrino  _et al._ [2022] F. M. D. Pellegrino, G. Falci, and E. Paladino, Effect of dilute impurities on short graphene josephson junctions, Communications Physics 5, 10.1038/s42005-022-01042-7 (2022).
  * Stosiek  _et al._ [2022] M. Stosiek, C. Baretzky, T. Balashov, F. Evers, and W. Wulfhekel, Friedel oscillations and superconducting gap enhancement by impurity scattering, Phys. Rev. B 105, L140504 (2022).
  * Sato and Ando [2017] M. Sato and Y. Ando, Topological superconductors: a review, Reports on Progress in Physics 80, 076501 (2017).
  * Feldman  _et al._ [2016] B. E. Feldman, M. T. Randeria, J. Li, S. Jeon, Y. Xie, Z. Wang, I. K. Drozdov, B. A. Bernevig, and A. Yazdani, High-resolution studies of the majorana atomic chain platform, Nature Physics 13, 286 (2016).
  * Schneider  _et al._ [2020] L. Schneider, S. Brinker, M. Steinbrecher, J. Hermenau, T. Posske, M. dos Santos Dias, S. Lounis, R. Wiesendanger, and J. Wiebe, Controlling in-gap end states by linking nonmagnetic atoms and artificially-constructed spin chains on superconductors, Nature Communications 11, 10.1038/s41467-020-18540-3 (2020).
  * Schneider  _et al._ [2021] L. Schneider, P. Beck, T. Posske, D. Crawford, E. Mascot, S. Rachel, R. Wiesendanger, and J. Wiebe, Topological shiba bands in artificial spin chains on superconductors, Nature Physics 17, 943 (2021).
  * Schneider  _et al._ [2022] L. Schneider, P. Beck, J. Neuhaus-Steinmetz, L. Rózsa, T. Posske, J. Wiebe, and R. Wiesendanger, Precursors of majorana modes and their length-dependent energy oscillations probed at both ends of atomic shiba chains, Nature Nanotechnology 17, 384 (2022).
  * Brinker  _et al._ [2022] S. Brinker, F. Küster, S. S. P. Parkin, P. Sessi, and S. Lounis, Anomalous excitations of atomically crafted quantum magnets, Science Advances 8, eabi7291 (2022), https://www.science.org/doi/pdf/10.1126/sciadv.abi7291 .
  * Küster  _et al._ [2022] F. Küster, S. Brinker, R. Hess, D. Loss, S. S. P. Parkin, J. Klinovaja, S. Lounis, and P. Sessi, Non-majorana modes in diluted spin chains proximitized to a superconductor, Proceedings of the National Academy of Sciences 119, 10.1073/pnas.2210589119 (2022).
  * Rüßmann and Blügel [2022] P. Rüßmann and S. Blügel, Density functional bogoliubov-de gennes analysis of superconducting nb and nb(110) surfaces, Phys. Rev. B 105, 125143 (2022).
  * Beck  _et al._ [2021] P. Beck, L. Schneider, L. Rózsa, K. Palotás, A. Lászlóffy, L. Szunyogh, J. Wiebe, and R. Wiesendanger, Spin-orbit coupling induced splitting of yu-shiba-rusinov states in antiferromagnetic dimers, Nature Communications 12, 10.1038/s41467-021-22261-6 (2021).
  * Mier  _et al._ [2021] C. Mier, D.-J. Choi, and N. Lorente, Calculations of in-gap states of ferromagnetic spin chains on $s$-wave wide-band superconductors, Phys. Rev. B 104, 245415 (2021).
  * Lászlóffy  _et al._ [2023] A. Lászlóffy, B. Nyári, G. Csire, L. Szunyogh, and B. Újfalussy, Topological superconductivity from first principles. ii. effects from manipulation of spin spirals: Topological fragmentation, braiding, and quasi-majorana bound states, Phys. Rev. B 108, 134513 (2023).
  * Nyári  _et al._ [2023] B. Nyári, A. Lászlóffy, G. Csire, L. Szunyogh, and B. Újfalussy, Topological superconductivity from first principles. i. shiba band structure and topological edge states of artificial spin chains, Phys. Rev. B 108, 134512 (2023).
  * Rüßmann  _et al._ [2023] P. Rüßmann, M. Bahari, S. Blügel, and B. Trauzettel, Interorbital cooper pairing at finite energies in rashba surface states, Phys. Rev. Res. 5, 043181 (2023).
  * Oliveira  _et al._ [1988] L. N. Oliveira, E. K. U. Gross, and W. Kohn, Density-functional theory for superconductors, Phys. Rev. Lett. 60, 2430 (1988).
  * Lüders  _et al._ [2005] M. Lüders, M. A. L. Marques, N. N. Lathiotakis, A. Floris, G. Profeta, L. Fast, A. Continenza, S. Massidda, and E. K. U. Gross, Ab initio theory of superconductivity. i. density functional formalism and approximate functionals, Phys. Rev. B 72, 024545 (2005).
  * Marques  _et al._ [2005] M. A. L. Marques, M. Lüders, N. N. Lathiotakis, G. Profeta, A. Floris, L. Fast, A. Continenza, E. K. U. Gross, and S. Massidda, Ab initio theory of superconductivity. ii. application to elemental metals, Phys. Rev. B 72, 024546 (2005).
  * Essenberger  _et al._ [2016] F. Essenberger, A. Sanna, P. Buczek, A. Ernst, L. Sandratskii, and E. K. U. Gross, Ab initio theory of iron-based superconductors, Phys. Rev. B 94, 014503 (2016).
  * Khaydukov  _et al._ [2018] Y. N. Khaydukov, A. S. Vasenko, E. A. Kravtsov, V. V. Progliado, V. D. Zhaketov, A. Csik, Y. V. Nikitenko, A. V. Petrenko, T. Keller, A. A. Golubov, M. Y. Kupriyanov, V. V. Ustinov, V. L. Aksenov, and B. Keimer, Magnetic and superconducting phase diagram of nb/gd/nb trilayers, Phys. Rev. B 97, 144511 (2018).
  * Pal  _et al._ [2022] B. Pal, A. Chakraborty, P. K. Sivakumar, M. Davydova, A. K. Gopi, A. K. Pandeya, J. A. Krieger, Y. Zhang, M. Date, S. Ju, N. Yuan, N. B. M. Schröter, L. Fu, and S. S. P. Parkin, Josephson diode effect from cooper pair momentum in a topological semimetal, Nature Physics 18, 1228 (2022).
  * Odobesko  _et al._ [2019] A. B. Odobesko, S. Haldar, S. Wilfert, J. Hagen, J. Jung, N. Schmidt, P. Sessi, M. Vogt, S. Heinze, and M. Bode, Preparation and electronic properties of clean superconducting Nb (110) surfaces, Phys. Rev. B 99, 115437 (2019).
  * Küster  _et al._ [2021] F. Küster, S. Brinker, S. Lounis, S. S. P. Parkin, and P. Sessi, Long range and highly tunable interaction between local spins coupled to a superconducting condensate, Nature Communications 12, 10.1038/s41467-021-26802-x (2021).
  * Küster  _et al._ [2021] F. Küster, A. M. Montero, F. S. M. Guimarães, S. Brinker, S. Lounis, S. S. P. Parkin, and P. Sessi, Correlating Josephson supercurrents and Shiba states in quantum spins unconventionally coupled to superconductors, Nature Communications 12, 10.1038/s41467-021-21347-5 (2021).
  * Odobesko  _et al._ [2020] A. Odobesko, D. Di Sante, A. Kowalski, S. Wilfert, F. Friedrich, R. Thomale, G. Sangiovanni, and M. Bode, Observation of tunable single-atom yu-shiba-rusinov states, Phys. Rev. B 102, 174504 (2020).
  * Friedrich  _et al._ [2021] F. Friedrich, R. Boshuis, M. Bode, and A. Odobesko, Coupling of yu-shiba-rusinov states in one-dimensional chains of fe atoms on nb(110), Phys. Rev. B 103, 235437 (2021).
  * Soldini  _et al._ [2023] M. O. Soldini, F. Küster, G. Wagner, S. Das, A. Aldarawsheh, R. Thomale, S. Lounis, S. S. P. Parkin, P. Sessi, and T. Neupert, Two-dimensional shiba lattices as a possible platform for crystalline topological superconductivity, Nature Physics 10.1038/s41567-023-02104-5 (2023).
  * Goedecke  _et al._ [2022] J. J. Goedecke, L. Schneider, Y. Ma, K. T. That, D. Wang, J. Wiebe, and R. Wiesendanger, Correlation of magnetism and disordered shiba bands in fe monolayer islands on nb(110), ACS Nano 16, 14066 (2022), pMID: 36001503, https://doi.org/10.1021/acsnano.2c03965 .
  * Beck  _et al._ [2022] P. Beck, L. Schneider, L. Bachmann, J. Wiebe, and R. Wiesendanger, Structural and superconducting properties of ultrathin ir films on nb(110), Phys. Rev. Mater. 6, 024801 (2022).
  * Lo Conte  _et al._ [2022] R. Lo Conte, M. Bazarnik, K. Palotás, L. Rózsa, L. Szunyogh, A. Kubetzka, K. von Bergmann, and R. Wiesendanger, Coexistence of antiferromagnetism and superconductivity in mn/nb(110), Phys. Rev. B 105, L100406 (2022).
  * Souza Mendes Guimarães  _et al._ [2023] F. Souza Mendes Guimarães, J. R. Suckert, H. Hamamera, S. Brinker, U. A. Aceves Rodríguez, J. Bouaziz, M. dos Santos Dias, and S. Lounis, TITAN, Zenodo 10.5281/zenodo.8099072 (2023).
  * Guimarães  _et al._ [2017] F. S. M. Guimarães, M. dos Santos Dias, J. Bouaziz, A. T. Costa, R. B. Muniz, and S. Lounis, Dynamical amplification of magnetoresistances and hall currents up to the thz regime, Scientific Reports 7, 3686 (2017).
  * Guimarães  _et al._ [2019] F. S. M. Guimarães, J. R. Suckert, J. Chico, J. Bouaziz, M. dos Santos Dias, and S. Lounis, Comparative study of methodologies to compute the intrinsic gilbert damping: interrelations, validity and physical consequences, Journal of Physics: Condensed Matter 31, 255802 (2019).
  * Guimarães  _et al._ [2020] F. S. M. Guimarães, J. Bouaziz, M. dos Santos Dias, and S. Lounis, Spin-orbit torques and their associated effective fields from gigahertz to terahertz, Communications Physics 3, 19 (2020).
  * Papaconstantopoulos [2015] D. Papaconstantopoulos,  _Handbook of the Band Structure of Elemental Solids: From Z = 1 To Z = 112_ , 2nd ed. (Springer US, 2015).
  * Guimarães  _et al._ [2015] F. S. M. Guimarães, S. Lounis, A. T. Costa, and R. B. Muniz, Dynamical current-induced ferromagnetic and antiferromagnetic resonances, Phys. Rev. B 92, 220410 (2015).
  * Hamamera  _et al._ [2022] H. Hamamera, F. S. M. Guimarães, M. dos Santos Dias, and S. Lounis, Polarisation-dependent single-pulse ultrafast optical switching of an elementary ferromagnet, Communications Physics 5, 16 (2022).
  * Hamamera  _et al._ [2023] H. Hamamera, F. S. M. Guimarães, M. dos Santos Dias, and S. Lounis, Ultrafast light-induced magnetization in non-magnetic films: from orbital and spin Hall phenomena to the inverse Faraday effect, arXiv: 2312.07888, 10.48550/arXiv.2312.07888 (2023).
  * Aceves Rodríguez  _et al._ [2023] U. A. Aceves Rodríguez, F. S. Mendes Guimarães, S. Brinker, and S. Lounis, Magnetic exchange interactions at the proximity of a superconductor, arXiv: 2306.02906, 10.48550/arXiv.2306.02906 (2023).
  * Tinkham [2004] M. Tinkham,  _Introduction to Superconductivity_ , Dover Books on Physics Series (Dover Publications, 2004).
  * De Gennes and Pincus [1999b] P. De Gennes and P. Pincus,  _Superconductivity Of Metals And Alloys_ , Advanced Books Classics (Avalon Publishing, 1999).
  * Ikushima and Mizusaki [1969] A. Ikushima and T. Mizusaki, Superconductivity in niobium and niobium-tantalum alloys, Journal of Physics and Chemistry of Solids 30, 873 (1969).
  * Veedu  _et al._ [2014] S. Veedu, Z. Popovic, and S. Satpathy, Theoretical model for rashba spin-orbit interaction in d electrons, Physical Review B 90 (2014).
  * Csire  _et al._ [2016] G. Csire, S. Schönecker, and B. Újfalussy, First-principles approach to thin superconducting slabs and heterostructures, Phys. Rev. B 94, 140502 (2016).
  * Klein and Papaconstantopoulos [1974] B. M. Klein and D. A. Papaconstantopoulos, Electron-phonon interaction and superconductivity in transition metals and transition-metal carbides, Phys. Rev. Lett. 32, 1193 (1974).
  * Saunderson  _et al._ [2020b] T. G. Saunderson, J. F. Annett, B. Újfalussy, G. Csire, and M. Gradhand, Gap anisotropy in multiband superconductors based on multiple scattering theory, Phys. Rev. B 101, 064510 (2020b).
  * Li  _et al._ [2016] J. Li, T. Neupert, Z. Wang, A. H. MacDonald, A. Yazdani, and B. A. Bernevig, Two-dimensional chiral topological superconductivity in shiba lattices, Nature Communications 7, 12297 (2016).
  * Powell [1964] M. Powell, An efficient method for finding the minimum of a function of several variables without calculating derivatives, The Computer Journal 7, 155 (1964).
  * Kodama  _et al._ [1983] J. Kodama, M. Itoh, and H. Hirai, Superconducting transition temperature versus thickness of Nb film on various substrates, Journal of Applied Physics 54, 4050 (1983), https://pubs.aip.org/aip/jap/article-pdf/54/7/4050/10571063/4050_1_online.pdf .
  * Krause and Thörnig [2018] D. Krause and P. Thörnig, JURECA: Modular supercomputer at Jülich Supercomputing Centre, Journal of large-scale research facilities JLSRF 4, 10.17815/jlsrf-4-121-1 (2018).


# Convergence in distribution of the Bernstein-Durrmeyer kernel and pointwise
convergence of a generalised operator for functions of bounded variation

Mohammed Taariq Mowzer<EMAIL_ADDRESS>Department of Mathematics and
Applied Mathematics, University of Cape Town

###### Abstract

We study the convergence of Bernstein type operators leading to two results.
The first: The kernel $K_{n}$ of the Bernstein-Durrmeyer operator at each
point $x\in(0,1)$ — that is $K_{n}(x,t)dt$ — once standardised converges to
the normal distribution. The second result computes the pointwise limit of a
generalised Bernstein-Durrmeyer operator applied to — possibly discontinuous —
functions $f$ of bounded variation.

###### keywords:

Positive linear operators , Bernstein type operators , Beta distribution

##  1 Introduction

######  Definition 1.1.

The Bernstein-Durrmeyer operator $M_{n}$ applied to an integrable function $f$
on $[0,1]$ is given by

$M_{n}(f)(x)=\int_{0}^{1}f(t)K_{n}(x,t)\,dt,\quad x\in[0,1],$

where $K_{n}$ is a Bernstein-Durrmeyer kernel given by

$K_{n}(x,t)=(n+1)\sum_{k=0}^{n}p_{n,k}(x)p_{n,k}(t),\quad x,t\in[0,1].$

The $p_{n,k}$ are the so called Bernstein polynomials defined as

$p_{n,k}(t)=\binom{n}{k}t^{k}(1-t)^{n-k}.$

Bernstein operators $B_{n}(f)(x)=\sum_{k=0}^{n}p_{n,k}(x)f(k/n)$ were
introduced by S.Bernstein [4] to give a constructive proof of the Weierstrass
approximation theorem. The Bernstein-Durrmeyer operator above was constructed
by Durrmeyer [7] and independently by Lupaş [10].

In his original proof, Bernstein used the binomial distribution and the weak
law of large numbers to constructively prove the Weierstrass approximation
theorem.

In this paper, we use the fact that the binomial distribution converges to the
normal distribution as the number of trials increases to show that the kernel
of the Bernstein-Durrmeyer operator tends to the normal distribution. To be
precise:

######  Theorem 1.2.

For $x\in(0,1)$ fixed

$\frac{1}{\sqrt{n}}K_{n}\left(x,\frac{t}{\sqrt{n}}+x\right)dt\text{ converges
in distribution to }N(0,2x(1-x))\text{ as }n\rightarrow\infty.$

Equivalently, if $X_{n}$ is a random variable with distribution $K_{n}(x,t)dt$
then

$\sqrt{n}(X_{n}-x)\text{ converges in distribution to }N(0,2x(1-x))\text{ as
}n\rightarrow\infty.$

The Bernstein-Durrmeyer operator can be generalised to arbitrary measure.
Here, we limit ourselves to arbitrary weights.

######  Definition 1.3.

The Bernstein-Durrmeyer operator with non-negative weight $w$ is given by

$M_{n,w}(f)(x)=\sum_{k=0}^{n}p_{n,k}(x)\frac{\int_{0}^{1}f(t)p_{n,k}(t)w(t)\,dt}{\int_{0}^{1}p_{n,k}(t)w(t)\,dt},$
(1)

where $w$ is chosen so that the terms are well defined, $x\in[0,1]$, and $f$ a
function on $[0,1]$ integrable with respect to $w(t)dt$.

In the particular case when $w$ is a Jacobi weight
$w(t)=t^{\alpha}(1-t)^{\beta}$ for $\alpha,\beta>-1$, the operator $M_{n,w}$
is well known. For example, Berens and Xu (e.g. [3]) investigated the
Bernstein-Durrmeyer operator with Jacobi weights, and Ditzian (e.g. [6])
studied the multidimensional Bernstein-Durrmeyer operator with Jacobi weights.
The case of arbitrary measure — replace $w(t)dt$ with $d\rho(t)$ in (1) — and
its multivariate analogue was, to our knowledge, for the first time studied by
Berdysheva and Jetter [2]. Berysheva [1] showed that the Bernstein-Durrmeyer
operator with arbitrary measure when applied to $f$ at $x$ converges to $f(x)$
if $f$ is essentially bounded, continuous at $x$ and the support of $\rho$
contains $x$.

Herzog and Hill proved in [9] that for $f$ defined on $\mathbb{Q}\cap[0,1]$,
$f$ bounded and having left and right limits $f(x-)$ and $f(x+)$, for the
classical Bernstein operator

$\lim_{n\rightarrow\infty}B_{n}(f)(x)=\frac{1}{2}f(x-)+\frac{1}{2}f(x+).$

Guo [8] showed when $f$ is of bounded variation that

$\lim_{n\rightarrow\infty}M_{n}(f)(x)=\frac{1}{2}f(x-)+\frac{1}{2}f(x+)$ (2)

holds for the Bernstein-Durrmeyer operator and gave an estimate of the rate of
convergence.

In this paper we find the limit of the Bernstein-Durrmeyer operator with
arbitrary weights when applied to a function $f$ of bounded variation. To be
precise, we show:

######  Theorem 1.4.

Let $w$ be a non-negative function of bounded variation on $[0,1]$. Let
$x\in(0,1)$ be fixed. Assume that either the left limit $w(x-)$ or the right
limit $w(x+)$ at $x$ is non-zero. For any function $f$ of bounded variation on
$[0,1]$,

$\lim_{n\rightarrow\infty}M_{n,w}(f)(x)=(1-\nu)f(x-)+\nu f(x+),$ (3)

where

$\nu=\int_{0}^{1}\frac{w(x+)(1-u)}{w(x-)u+w(x+)(1-u)}\,du.$

When $w(x-)=w(x+)$, $\nu$ evaluates to $\frac{1}{2}$. Otherwise,

$\nu=\frac{r^{2}-r[1+\ln r]}{(r-1)^{2}},\quad\text{for
}r=\frac{w(x+)}{w(x-)}\text{ when }r\in(0,\infty)\text{ and }r\neq 1.$

Theorem 1.4 is an extension of (2) — consider $w=1$. In fact (2) holds even if
only $w(x+)=w(x-)$, that is, $w$ continuous at $x$.

##  2 Prerequisites and useful results

Henceforth, we denote the binomial distribution with $n$ trials and
probability $p$ of success by $Bin(n,p)$, the normal distribution with mean
$\mu$ and variance $\sigma^{2}$ by $N(\mu,\sigma^{2})$, and the beta
distribution with $\alpha$ and $\beta$ parameters by $Beta(\alpha,\beta)$. The
symbol $\stackrel{{\scriptstyle D}}{{\rightarrow}}$ denotes convergence in
distribution and $\sim$ denotes having the same distribution.

The function $p_{n,k}$ is extended onto $\mathbb{R}$, and $k$ is allowed on
$[0,n]$ so that

$p_{n,k}(t)=\begin{cases}\frac{\Gamma(n+1)}{\Gamma(k+1)\Gamma(n-k+1)}t^{k}(1-t)^{n-k}&\text{
if }t\in[0,1],\\\ 0&\text{ if }t\not\in[0,1].\end{cases}$

For a function $f$, $f(t-)$ and $f(t+)$ denote the left and right limit of $f$
at $t$, and $\mathcal{D}(f)$ is the set of discontinuities of $f$. Usually $x$
will be fixed in $(0,1)$. It is safe to define $\Psi$ as the cumulative
distribution function (cdf) of $N(0,x(1-x))$.

######  Theorem 2.1.

Let $C_{n}$ have a $Bin(n,p)$ distribution for fixed $p\in(0,1)$. Then

$\sqrt{n}(C_{n}/n-p)\stackrel{{\scriptstyle
D}}{{\rightarrow}}N(0,p(1-p)),\text{ as }n\rightarrow\infty.$

######  Theorem 2.2.

Let $B_{r_{1},r_{2}}$ have distribution $Beta(r_{1},r_{2})$. Then

$\sqrt{r_{1}+r_{2}}\left(B_{r_{1},r_{2}}-\frac{r_{1}}{r_{1}+r_{2}}\right)\stackrel{{\scriptstyle
D}}{{\rightarrow}}N(0,\gamma(1-\gamma)),$ (4)

as $r_{1},r_{2}\rightarrow\infty$ and
$\frac{r_{1}}{r_{1}+r_{2}}\rightarrow\gamma$ for $\gamma\in(0,1)$.

Theorem 2.1 is a direct consequence of the central limit theorem and is well
known. Theorem 2.2 is a generalisation of the well known approximation that
the standardised $Beta(an,bn)$ distribution tends to the normal distribution
as $n\rightarrow\infty$. As far as we know, Theorem 2.2 is not well known. We
supply a proof in the appendix.

The Corollary below will be the most useful form of Theorem 2.2.

######  Corollary 2.3.

Suppose $B_{n,k}$ has a $Beta(k+1,n-k+1)$ distribution, and let $x\in(0,1)$ be
fixed. Then

$\sqrt{n}\left(B_{n,k}-x\right)\stackrel{{\scriptstyle
D}}{{\rightarrow}}N(\alpha,x(1-x))$ (5)

as $n\rightarrow\infty$ and
$\sqrt{n}\left(\frac{k}{n}-x\right)\rightarrow\alpha$.

###### Proof.

Notice $\sqrt{n}\left(\frac{k}{n}-x\right)\rightarrow\alpha$ implies
$\frac{k}{n}\rightarrow x$. By Theorem 2.2,

${\overbrace{\frac{\sqrt{n+2}}{\sqrt{n}}}^{M}}\left[\sqrt{n}\left(B_{n,k}-\frac{k}{n}\right)\right]+{\overbrace{\sqrt{n+2}\left(\frac{k}{n}-\frac{k+1}{n+2}\right)}^{C}}=\sqrt{n+2}\left(B_{n,k}-\frac{k+1}{n+2}\right)\stackrel{{\scriptstyle
D}}{{\rightarrow}}N(0,x(1-x))$

as $n\rightarrow\infty$ and $\frac{k}{n}\rightarrow x$ . Notice $M\rightarrow
1$ and $C\rightarrow 0$ in the limit. Therefore, by Slutsky’s Theorem,
$\sqrt{n}\left(B_{n,k}-\frac{k}{n}\right)\stackrel{{\scriptstyle
D}}{{\rightarrow}}N(0,x(1-x))$. Another application of Slutsky’s Theorem gives

$\sqrt{n}(B_{n,k}-x)=\sqrt{n}\left(B_{n,k}-\frac{k}{n}\right)+\sqrt{n}\left(\frac{k}{n}-x\right)\stackrel{{\scriptstyle
D}}{{\rightarrow}}N(\alpha,x(1-x))$

as $n\rightarrow\infty$ and
$\sqrt{n}\left(\frac{k}{n}-x\right)\rightarrow\alpha$. ∎

The following gives a straightforward way to show uniform convergence on
compact intervals using sequences.

######  Lemma 2.4.

[11, Chapter 3]. Let $f$ and $(f_{n})$ be functions defined on a compact
subset $K$ of $\mathbb{R}$. The following are equivalent

  1. i)

$f_{n}$ converges uniformly to $f$, and $f$ is continuous on $K$.

  2. ii)

For any sequence $(\alpha_{n})\subset K$, if
$\lim_{n\rightarrow\infty}a_{n}=a$ then
$\lim_{n\rightarrow\infty}f_{n}(\alpha_{n})=f(\alpha)$.

###### Proof.

See for example [11, Chapter 3] section on continuous convergence. ∎

The next Theorem allows us to use convergence in distribution to show
convergence in expectation.

######  Theorem 2.5.

Let $(X_{n})$ be a sequence of real-valued random variables which converge to
the random variable $X$ in distribution. Let $(h_{n})$ be a sequence of
uniformly bounded measurable functions from $\mathbb{R}$ to $\mathbb{R}$ which
converge uniformly on compact intervals to a measurable function $h$ whose set
of discontinuities $\mathcal{D}(h)$ is closed and satisfies
$\mathbb{P}(X\not\in\mathcal{D}(h))=1$. Then

$\lim_{n\rightarrow\infty}\mathbb{E}(h_{n}(X_{n}))=\mathbb{E}(h(X)).$

###### Proof.

By Skorokhod’s representation theorem, see for example [5, Chapter 3.5], there
exists a sequence of random variables $(Y_{n})$ such that each $Y_{n}$ has the
same distribution as $X_{n}$, and the $Y_{n}$ converges almost surely to a
random variable $Y$ having the same distribution as $X$. Let $\Omega_{0}$ be
the almost sure set where $(Y_{n})$ converges pointwise to $Y$ and let
$\Omega_{1}$ be the almost sure set $\\{Y\not\in\mathcal{D}(h)\\}$. For
$\omega$ in the almost sure set $\Omega_{0}\cap\Omega_{1}$, $Y(\omega)$ is in
the open set $\mathcal{D}(h)^{c}$. Hence, $h$ is continuous in a compact
neighbourhood of $Y(\omega)$. By Lemma 2.4, uniform convergence of $h_{n}$ to
$h$ on the compact neighbourhood yields $h_{n}(Y_{n}(\omega))\rightarrow
h(Y(\omega))$ as $n\rightarrow\infty$. Therefore, $h_{n}(Y_{n})\rightarrow
h(Y)$ almost surely. Since $(h_{n})$ is uniformly bounded, the sequence
$(h_{n}(Y_{n}))$ is dominated by some constant. By Lebesgue’s dominated
convergence theorem,

$\lim_{n\rightarrow\infty}\mathbb{E}(h_{n}(X_{n}))=\lim_{n\rightarrow\infty}\mathbb{E}(h_{n}(Y_{n}))=\mathbb{E}(h(Y))=\mathbb{E}(h(X)).$

∎

##  3 Normality of Bernstein-Durrmeyer kernel

For the rest of this section $x\in(0,1)$ is fixed. We construct a
probabilistic representation of the Bernstein-Durrmeyer kernel
$K_{n}(x,\cdot)$. Consider the distribution generated from the samples of the
following process.

  1. Step 1.

Obtain a coin with probability $x$ of obtaining heads.

  2. Step 2.

Flip the coin $n$ times recording the number of heads.

  3. Step 3.

If $k$ heads were recorded, the sample is randomly drawn from a
$Beta(k+1,n-k+1)$ distribution.

Suppose $C_{n}\sim Bin(n,x)$ and
$B_{n,k}~{}\sim~{}Beta(k~{}+~{}1,~{}n~{}-~{}k~{}+~{}1)$ form an independency
of random variables. The random variable

$X_{n}=\sum_{k=0}^{n}\mbox{\large$\chi$}_{\\{C_{n}=k\\}}B_{n,k}$

describes samples from the above process. Here, $\mbox{\large$\chi$}_{F}$ is
the indicator function on a set $F$ that is

$\mbox{\large$\chi$}_{F}(x)=\begin{cases}1&\text{ if }x\in F,\\\ 0&\text{ if
}x\not\in F.\end{cases}$

We show $X_{n}$ has distribution $K_{n}(x,\cdot)dt$. First, we note the
probability density function (pdf) of $B_{n,k}$ is

$f_{B_{n,k}}(t)=\frac{\Gamma(n+2)}{\Gamma(k+1)\Gamma(n-k+1)}t^{k}(1-t)^{n-k}=(n+1)p_{n,k}(t)\quad\text{for
}t\in[0,1]$

and the probability mass function (pmf) of $C_{n}$ is

$p_{C_{n}}(k)=\binom{n}{k}x^{k}(1-x)^{n-k}=p_{n,k}(x)\quad\text{ for
}k\in\\{0,\ldots,n\\}.$

Now, for all $z\in\mathbb{R}$,

$\displaystyle\mathbb{P}(X_{n}\leq z)$
$\displaystyle=\mathbb{E}[\mathbb{P}(X_{n}\leq z\,|\,C_{n}=k)]$
$\displaystyle=\mathbb{E}[\mathbb{P}(B_{n,k}\leq z\,|\,C_{n}=k)]$
$\displaystyle=\sum_{k=0}^{n}p_{n,k}(x)\mathbb{P}(B_{n,k}\leq z)$
$\displaystyle=\sum_{k=0}^{n}p_{n,k}(x)\int_{0}^{z}(n+1)p_{n,k}(t)\,dt$
$\displaystyle=\int_{0}^{z}K_{n}(x,t)\,dt.$

Independence of $C_{n}$ and $B_{n,k}$ is used in the third equality.

######  Theorem 3.1.

Let $x\in(0,1)$ be fixed. Suppose $X_{n}$ is defined as above. Then

$\sqrt{n}(X_{n}-x)\stackrel{{\scriptstyle D}}{{\rightarrow}}N(0,2x(1-x)).$

###### Proof.

Let $z\in\mathbb{R}$ be fixed. Let $C_{n}^{*}=\sqrt{n}(C_{n}/n-x)$ and

$h_{n}(\alpha)=\begin{cases}\mathbb{P}(\sqrt{n}(B_{n,xn+\alpha\sqrt{n}}-x)\leq
z)&\text{ if }xn+\alpha\sqrt{n}\in[0,n],\\\ 0&\text{ otherwise.}\end{cases}$

Then

$\displaystyle\mathbb{P}(\sqrt{n}(X_{n}-x)\leq z)$
$\displaystyle=\mathbb{E}[\mathbb{P}(\sqrt{n}(B_{n,k}-x)\leq z\,|\,C_{n}=k)]$
$\displaystyle=\mathbb{E}[\mathbb{P}(\sqrt{n}(B_{n,xn+\alpha\sqrt{n}}-x)\leq
z\,|\,\sqrt{n}(C_{n}/n-x)=\alpha)]$
$\displaystyle=\mathbb{E}[\mathbb{P}(\sqrt{n}(B_{n,xn+\alpha\sqrt{n}}-x)\leq
z\,|\,C_{n}^{*}=\alpha)]=\mathbb{E}(h_{n}(C_{n}^{*})),$ (6)

where the last equality follows from the independence of $C_{n}$ and
$B_{n,k}$.

Let $K$ be any non-empty compact subset of $\mathbb{R}$. Let
$(\alpha_{n})\subset K$ such that
$\lim_{n\rightarrow\infty}\alpha_{n}=\alpha$. Then
$\lim_{n\rightarrow\infty}\sqrt{n}\left(\frac{xn+\alpha_{n}\sqrt{n}}{n}-x\right)=\alpha$.
Corollary 2.3 then gives
$\sqrt{n}(B_{n,xn+\alpha_{n}\sqrt{n}}-x)\stackrel{{\scriptstyle
D}}{{\rightarrow}}N(\alpha,x(1-x))$ as $n\rightarrow\infty$. Let $A$ be a
random variable with distribution $N(\alpha,x(1-x))$ and recall $\Psi$ is the
cdf of $N(0,x(1-x))$. Convergence in distribution yields,

$\lim_{n\rightarrow\infty}h_{n}(\alpha_{n})=\lim_{n\rightarrow\infty}\mathbb{P}(\sqrt{n}(B_{n,xn+\alpha_{n}\sqrt{n}}-x)\leq
z)=\mathbb{P}(A\leq z)=\Psi(z-\alpha).$ (7)

Lemma 2.4 gives that $h_{n}$ converges uniformly to $\Psi(\cdot-\alpha)$ on
$K$. Since $K$ was arbitrary, uniform convergence holds on all compact
intervals.

Let $U$ and $V$ be independent random variables with distribution
$N(0,x(1-x)).$ By Theorem 2.1, $C_{n}^{*}\stackrel{{\scriptstyle
D}}{{\rightarrow}}V$ as $n\rightarrow\infty$ . Equations (6), (7) and Theorem
2.5 give

$\displaystyle\lim_{n\rightarrow\infty}\mathbb{P}(\sqrt{n}(X_{n}-x)\leq
z)=\lim_{n\rightarrow\infty}\mathbb{E}(h_{n}(C_{n}^{*}))=\mathbb{E}(\Psi(z-V))$
$\displaystyle=\mathbb{E}(\mathbb{P}(U\leq z-V\,|\,V))$
$\displaystyle=\mathbb{P}(U\leq z-V)=\mathbb{P}(U+V\leq z).$

Since $U$ and $V$ are independent, the sum $U+V$ has distribution
$N(0,2x(1-x))$. Hence, the cdf of $\sqrt{n}(X_{n}-x)$ converges pointwise to
the cdf of $N(0,2x(1-x))$. ∎

##  4 Limit of the generalised Bernstein-Durrmeyer operator

For the rest of this section $x\in(0,1)$ is fixed. Recall,

$M_{n,w}(f)(x)=\sum_{k=0}^{n}p_{n,k}(x)\frac{\int_{0}^{1}f(t)p_{n,k}(t)w(t)\,dt}{\int_{0}^{1}p_{n,k}(t)w(t)\,dt}.$

The following lemma evaluates the limit of the integral appearing in the
equation above.

######  Lemma 4.1.

Let $f:[0,1]\rightarrow\mathbb{R}$ be of bounded variation. Let
$k_{n}(\alpha)=xn+\alpha\sqrt{n}$. The limit

$\lim_{n\rightarrow\infty}(n+1)\int_{0}^{1}f(t)p_{n,k_{n}(\alpha)}(t)\,dt=\Psi(-\alpha)f(x-)+[1-\Psi(-\alpha)]f(x+),$
(8)

holds uniformly in $\alpha$ on all compact subsets of $\mathbb{R}$.

###### Proof.

We look to apply Lemma 2.4. Let $K$ be a non-empty compact subset of
$\mathbb{R}$. Let $(\alpha_{n})\subset K$ such that
$\lim_{n\rightarrow\infty}\alpha_{n}=\alpha$. Clearly,
$\lim_{n\rightarrow\infty}\sqrt{n}\left(\frac{k_{n}(\alpha_{n})}{n}-x\right)=\alpha$.
Hence, Corollary 2.3 yields

$B_{n,k_{n}(\alpha_{n})}^{*}:=\sqrt{n}\left(B_{n,k_{n}(\alpha_{n})}-x\right)\stackrel{{\scriptstyle
D}}{{\rightarrow}}N(\alpha,x(1-x))\text{ as }n\rightarrow\infty.$ (9)

Notice

$(n+1)\int_{0}^{1}f(t)p_{n,k_{n}(\alpha_{n})}(t)\,dt=\mathbb{E}(f(B_{n,k_{n}(\alpha_{n})}))=\mathbb{E}(h_{n}(B_{n,k_{n}(\alpha_{n})}^{*})),$
(10)

where

$h_{n}(t)=f\left(\frac{t}{\sqrt{n}}+x\right)\mbox{\large$\chi$}_{[0,1]}\left(\frac{t}{\sqrt{n}}+x\right)=f\left(\frac{t}{\sqrt{n}}+x\right)\mbox{\large$\chi$}_{[-x\sqrt{n},(1-x)\sqrt{n}]}(t),\quad
t\in\mathbb{R}.$

From the decomposition

$\mbox{\large$\chi$}_{[-x\sqrt{n},(1-x)\sqrt{n}]}=\mbox{\large$\chi$}_{[-x\sqrt{n},0)}+\mbox{\large$\chi$}_{\\{0\\}}+\mbox{\large$\chi$}_{(0,(1-x)\sqrt{n}]},$

we notice $h_{n}$ limits uniformly to
$h=f(x-)\mbox{\large$\chi$}_{(-\infty,0)}+f(x)\mbox{\large$\chi$}_{\\{0\\}}+f(x+)\mbox{\large$\chi$}_{(0,\infty)}$
on $K$. Since $f$ is of bounded variation, $(h_{n})$ is uniformly bounded. Let
$A$ be a random variable with distribution $N(\alpha,x(1-x))$. The set of
discontinuities $\mathcal{D}(h)\subset\\{0\\}$ is closed, and is a null set of
the distribution of $A$. Hence, equations (10), (9) and Theorem 2.5 give

$(n+1)\int_{0}^{1}f(t)p_{n,k_{n}(\alpha_{n})}(t)\,dt=\mathbb{E}(h_{n}(B_{n,k_{n}(\alpha_{n})}^{*}))\stackrel{{\scriptstyle
n\rightarrow\infty}}{{\longrightarrow}}\mathbb{E}(h(A))=\Psi(-\alpha)f(x-)+[1-\Psi(-\alpha)]f(x+).$
(11)

Lemma 2.4 applied to (11) says that convergence of (8) is uniform in $\alpha$
on $K$. Since $K$ was arbitrary (8) holds uniformly in $\alpha$ on all compact
subsets of $\mathbb{R}$. ∎

###### Proof of Theorem 1.4.

Note,

$M_{n,w}(f)(x)=\sum_{k=0}^{n}p_{n,k}(x)c_{n}(k)=\mathbb{E}(c_{n}(C_{n})),$
(12)

where $C_{n}\sim Bin(n,x)$ and

$c_{n}(k)=\frac{\int_{0}^{1}f(t)p_{n,k}(t)w(t)\,dt}{\int_{0}^{1}p_{n,k}(t)w(t)\,dt}\text{
if }k\in[0,n],\text{ otherwise }c_{n}(k)=0.$

We transform $C_{n}$ and $c_{n}$ as follows. Let
$C_{n}^{*}=\sqrt{n}\left(\frac{C_{n}}{n}-x\right)$ and define $c_{n}^{*}$ so
that $c_{n}(C_{n})=c_{n}^{*}(C_{n}^{*})$, that is
$c_{n}^{*}(\alpha)=c_{n}(xn+\alpha\sqrt{n})$ for $\alpha\in\mathbb{R}$. Then

$M_{n,w}(f)(x)=\mathbb{E}(c_{n}(C_{n}))=\mathbb{E}(c_{n}^{*}(C_{n}^{*})).$
(13)

We define functions $t_{n},b_{n}$ as follows

$t_{n}(\alpha)=(n+1)\int_{0}^{1}f(t)p_{n,xn+\alpha\sqrt{n}}(t)w(t)\,dt\quad\text{and}\quad
b_{n}(\alpha)=(n+1)\int_{0}^{1}p_{n,xn+\alpha\sqrt{n}}(t)w(t)\,dt$

for $\alpha\in\mathbb{R}$.

By Lemma 4.1,

$\lim_{n\rightarrow\infty}t_{n}(\alpha)=f(x-)w(x-)\Psi(-\alpha)+f(x+)w(x+)(1-\Psi(-\alpha))$
(14)

and

$\lim_{n\rightarrow\infty}b_{n}(\alpha)=w(x-)\Psi(-\alpha)+w(x+)(1-\Psi(-\alpha))$
(15)

uniformly in $\alpha$ on compact subsets of $\mathbb{R}$. By (14) and (15)

$\lim_{n\rightarrow\infty}c_{n}^{*}(\alpha)=\lim_{n\rightarrow\infty}\frac{t_{n}(\alpha)}{b_{n}(\alpha)}=\frac{f(x-)w(x-)\Psi(-\alpha)+f(x+)w(x+)(1-\Psi(-\alpha))}{w(x-)\Psi(-\alpha)+w(x+)(1-\Psi(-\alpha))}$
(16)

uniformly on compact subsets of $\mathbb{R}$. Notice $w(x-)\neq 0$ or
$w(x+)\neq 0$ is used here to guarantee that for every compact subset of
$\mathbb{R}$ there exists some $m>0$ which is a lower bound of the
denominator.

Since $f$ is of bounded variation it is bounded, so $(c_{n})$ and
$(c_{n}^{*})$ are uniformly bounded. Let $V$ have a $N(0,x(1-x))$
distribution. By Theorem 2.1, $(C_{n}^{*})$ converges in distribution to $V$.
Equation (16) says that $(c_{n}^{*})$ converges uniformly on compact intervals
to the right hand side of (16) which we shall call $c^{*}$. Recalling (13), an
application of Theorem 2.5 gives

$\lim_{n\rightarrow\infty}M_{n,w}(f)(x)=\lim_{n\rightarrow\infty}\mathbb{E}(c_{n}^{*}(C_{n}^{*}))=\mathbb{E}(c^{*}(V)).$
(17)

Expanding out the right hand side gives

$\lim_{n\rightarrow\infty}M_{n,w}(f)(x)=\int_{-\infty}^{\infty}\frac{f(x-)w(x-)\Psi(-\alpha)+f(x+)w(x+)(1-\Psi(-\alpha))}{w(x-)\Psi(-\alpha)+w(x+)(1-\Psi(-\alpha))}\,d\Psi(\alpha),$
(18)

The substitution $u=\Psi(-\alpha)$ in (18) gives (3). ∎

##  5 Acknowledgements

Much thanks to Professor E.E. Berdysheva for suggesting this problem, for
greatly helping in the preparation of the paper, and for providing reading
materials on the history of the Bernstein-Durrmeyer operators. The author also
thanks F.C. Rath for bringing to attention the work of Professor R.J. Ryder
[12] who the author thanks for providing further references on the normal
approximation of the beta distribution.

## References

  * [1] Berdysheva, E. E. Bernstein–Durrmeyer operators with respect to arbitrary measure, II: Pointwise convergence.  Journal of Mathematical Analysis and Applications 418, 2 (2014), 734–752.
  * [2] Berdysheva, E. E., and Jetter, K. Multivariate Bernstein–Durrmeyer operators with arbitrary weight functions.  Journal of Approximation Theory 162, 3 (2010), 576–598.  Special Issue: Bommerholz Proceedings.
  * [3] Berens, H., and Xu, Y. On Bernstein-Durrmeyer polynomials with Jacobi weights.  C. K. Chui, Ed., Academic Press Boston, pp. 25–46.
  * [4] Bernstein, S. Démonstration du théorème de weierstrass fondée sur le calcul des probabilités.  Communications de la Société mathématique de Kharkow. 2-ée série 13 (1912), 1–2.
  * [5] Çınlar, E. Probability and Stochastics.  Graduate Texts in Mathematics. Springer New York, 2011.
  * [6] Ditzian, Z. Multidimensional Jacobi-type Bernstein-Durrmeyer operators.  Acta Scientiarum Mathematicarum 60, 1 (1995), 225–244.
  * [7] Durrmeyer, J. L. Une formule d’inversion de la transformée de Laplace: Applications à la théorie des moments.  PhD thesis, Faculté des Sciences de l’Université de Paris, 1967\.
  * [8] Guo, S. On the rate of convergence of the Durrmeyer operator for functions of bounded variation.  Journal of Approximation Theory 51, 2 (1987), 183–192.
  * [9] Herzog, F., and Hill, J. D. The Bernstein polynomials for discontinuous functions.  American Journal of Mathematics 68, 1 (1946), 109–124.
  * [10] Lupaş, A. Die Folge der Betaoperatoren.  PhD thesis, Universität Stuttgart, 1972.
  * [11] Remmert, R., and Burckel, R. Theory of Complex Functions.  Graduate Texts in Mathematics. Springer New York, 1991.
  * [12] Ryder, R. J. The Beta$(b,b)$ distribution converges to the normal distribution when $b\rightarrow\infty$.  Univariate Distribution Relationships, 2013.
  * [13] Scheffe, H. A useful convergence theorem for probability distributions.  The Annals of Mathematical Statistics 18, 3 (1947), 434–438.

##  Appendix A Proof of Theorem 2.2

###### Proof.

We follow the proof method used in [12] which showed that the $Beta(n,n)$ goes
to the normal distribution as $n\rightarrow\infty$. That is for
$B_{r_{1},r_{2}}\sim Beta(r_{1},r_{2})$, we show that the pdf of
$\sqrt{r_{1}+r_{2}}(B_{r_{1},r_{2}}-\frac{r_{1}}{r_{1}+r_{2}})$ converges
pointwise to the pdf of the $N(0,\gamma(1-\gamma))$ distribution as
$r_{1},r_{2}\rightarrow\infty$ and
$\frac{r_{1}}{r_{1}+r_{2}}\rightarrow\gamma$. By Scheffe’s Theorem [13],
pointwise convergence of the pdf implies convergence in distribution.

The beta distribution has pdf

$f_{B_{r_{1},r_{2}}}(x)=\frac{\Gamma(r_{1}+r_{2})}{\Gamma(r_{1})\Gamma(r_{2})}x^{r_{1}-1}(1-x)^{r_{2}-1}\quad\text{for
}0<x<1.$

Using the transformation formula, the pdf of
$Y_{r_{1},r_{2}}=\sqrt{r_{1}+r_{2}}(B_{r_{1},r_{2}}-\frac{r_{1}}{r_{1}+r_{2}})$
is

$\displaystyle f_{Y_{r_{1},r_{2}}}(y)$
$\displaystyle=\frac{1}{\sqrt{r_{1}+r_{2}}}f_{B_{r_{1},r_{2}}}\left(\frac{r_{1}}{r_{1}+r_{2}}+\frac{y}{\sqrt{r_{1}+r_{2}}}\right)$
$\displaystyle=\frac{1}{\sqrt{r_{1}+r_{2}}}\frac{\Gamma(r_{1}+r_{2})}{\Gamma(r_{1})\Gamma(r_{2})}\left(\frac{r_{1}}{r_{1}+r_{2}}+\frac{y}{\sqrt{r_{1}+r_{2}}}\right)^{r_{1}-1}\left(\frac{r_{2}}{r_{1}+r_{2}}-\frac{y}{\sqrt{r_{1}+r_{2}}}\right)^{r_{2}-1},$

for $-\frac{r_{1}}{\sqrt{r_{1}+r_{2}}}<y<\frac{r_{2}}{\sqrt{r_{1}+r_{2}}}$. We
show pointwise convergence in the log of the pdf. Note

$\displaystyle\ln(f_{Y_{r_{1},r_{2}}}(y))=$
$\displaystyle\underbrace{-\ln{\sqrt{r_{1}+r_{2}}}+\ln(\Gamma(r_{1}+r_{2}))-\ln(\Gamma(r_{1}))-\ln(\Gamma(r_{2}))}_{A}$
$\displaystyle+\underbrace{(r_{1}-1)\ln\left(\frac{r_{1}}{r_{1}+r_{2}}+\frac{y}{\sqrt{r_{1}+r_{2}}}\right)+(r_{2}-1)\ln\left({\frac{r_{2}}{r_{1}+r_{2}}-\frac{y}{\sqrt{r_{1}+r_{2}}}}\right)}_{B}$

Stirling’s approximation gives $\ln(\Gamma(n))=n\ln
n-n+\ln\sqrt{2\pi/n}+O(\frac{1}{n})$, $n\rightarrow\infty$. Approximating $A$,

$\displaystyle A$
$\displaystyle=-\ln\sqrt{r_{1}+r_{2}}+(r_{1}+r_{2})\ln(r_{1}+r_{2})-(r_{1}+r_{2})+\ln\sqrt{2\pi/(r_{1}+r_{2})}+O\left(\frac{1}{r_{1}+r_{2}}\right)$
$\displaystyle\quad-\left[r_{1}\ln
r_{1}-r_{1}+\ln\sqrt{2\pi/r_{1}}+O\left(\frac{1}{r_{1}}\right)\right]-\left[r_{2}\ln
r_{2}-r_{2}+\ln\sqrt{2\pi/r_{2}}+O\left(\frac{1}{r_{2}}\right)\right]$
$\displaystyle=-\ln\sqrt{2\pi}+r_{1}\ln\frac{r_{1}+r_{2}}{r_{1}}+r_{2}\ln\frac{r_{1}+r_{2}}{r_{2}}+\frac{1}{2}\ln\frac{r_{1}}{r_{1}+r_{2}}+\frac{1}{2}\ln\frac{r_{2}}{r_{1}+r_{2}}+O\left(\frac{1}{r_{1}}+\frac{1}{r_{2}}\right).$

The limit is pointwise; so fix $y$. Note
$\frac{r_{1}}{r_{1}+r_{2}}\rightarrow\gamma$, so
$y\frac{\sqrt{r_{1}+r_{2}}}{r_{1}}$ and $y\frac{\sqrt{r_{1}+r_{2}}}{r_{2}}$
limit to 0. Hence, we can use the second order Taylor Series approximation
$\ln(1+x)=x-\frac{x^{2}}{2}+O(x^{3})$ as $x\rightarrow 0$. Approximating $B$,

$\displaystyle B$
$\displaystyle=(r_{1}-1)\left[\ln\frac{r_{1}}{r_{1}+r_{2}}+\ln(1+y\frac{\sqrt{r_{1}+r_{2}}}{r_{1}})\right]$
$\displaystyle\quad+(r_{2}-1)\left[\ln\frac{r_{2}}{r_{1}+r_{2}}+\ln(1-y\frac{\sqrt{r_{1}+r_{2}}}{r_{2}})\right]$
$\displaystyle=r_{1}\ln\frac{r_{1}}{r_{1}+r_{2}}+y\sqrt{r_{1}+r_{2}}-\frac{1}{2}y^{2}\frac{r_{1}+r_{2}}{r_{1}}-\ln\frac{r_{1}}{r_{1}+r_{2}}+O\left(\frac{1}{\sqrt{r_{1}}}\right)$
$\displaystyle\quad+r_{2}\ln\frac{r_{2}}{r_{1}+r_{2}}-y\sqrt{r_{1}+r_{2}}-\frac{1}{2}y^{2}\frac{r_{1}+r_{2}}{r_{2}}-\ln\frac{r_{2}}{r_{1}+r_{2}}+O\left(\frac{1}{\sqrt{r_{2}}}\right).$

Combining approximations for $A$ and $B$ we get

$\displaystyle A+B$
$\displaystyle=\left[r_{1}\ln\frac{r_{1}+r_{2}}{r_{1}}+r_{1}\ln\frac{r_{1}}{r_{1}+r_{2}}\right]+\left[r_{2}\ln\frac{r_{1}+r_{2}}{r_{2}}+r_{2}\ln\frac{r_{2}}{r_{1}+r_{2}}\right]+\left[y\sqrt{r_{1}+r_{2}}-y\sqrt{r_{1}+r_{2}}\right]$
$\displaystyle\quad+\left[\frac{1}{2}\ln\frac{r_{1}}{r_{1}+r_{2}}-\ln\frac{r_{1}}{r_{1}+r_{2}}\right]+\left[\frac{1}{2}\ln\frac{r_{2}}{r_{1}+r_{2}}-\ln\frac{r_{2}}{r_{1}+r_{2}}\right]$
$\displaystyle\quad-\ln\sqrt{2\pi}-\frac{1}{2}y^{2}\left(\frac{r_{1}+r_{2}}{r_{1}}+\frac{r_{1}+r_{2}}{r_{2}}\right)+O\left(\frac{1}{\sqrt{r_{1}}}+\frac{1}{\sqrt{r_{2}}}\right)$
$\displaystyle=\quad-\frac{1}{2}\ln\frac{r_{1}}{r_{1}+r_{2}}-\frac{1}{2}\ln\frac{r_{2}}{r_{1}+r_{2}}-\ln\sqrt{2\pi}-\frac{1}{2}y^{2}\left(\frac{r_{1}+r_{2}}{r_{1}}+\frac{r_{1}+r_{2}}{r_{2}}\right)+O\left(\frac{1}{\sqrt{r_{1}}}+\frac{1}{\sqrt{r_{2}}}\right)$
$\displaystyle\rightarrow-\frac{1}{2}\ln\gamma-\frac{1}{2}\ln(1-\gamma)-\ln\sqrt{2\pi}-\frac{1}{2}\frac{y^{2}}{\gamma(1-\gamma)}=\ln\left(\frac{1}{\sqrt{2\pi\gamma(1-\gamma)}}e^{-\frac{1}{2}\frac{y^{2}}{\gamma(1-\gamma)}}\right),$

as $r_{1},r_{2}\rightarrow\infty$ and
$\frac{r_{1}}{r_{1}+r_{2}}\rightarrow\gamma$. This is the natural log of the
pdf of the $N(0,\gamma(1-\gamma))$ distribution. ∎