# Blind Knowledge Distillation for Robust Image Classification

Timo Kaiser  Lukas Ehmann  Christoph Reinders&Bodo Rosenhahn Institute for
Information Processing, Leibniz University Hannover {kaiser, ehmannlu,
reinders<EMAIL_ADDRESS>
###### Abstract

Optimizing neural networks with noisy labels is a challenging task, especially
if the label set contains real-world noise. Networks tend to generalize to
reasonable patterns in the early training stages and overfit to specific
details of noisy samples in the latter ones. We introduce Blind Knowledge
Distillation \- a novel teacher-student approach for learning with noisy
labels by masking the ground truth related teacher output to filter out
potentially corrupted ‘knowledge’ and to estimate the tipping point from
generalizing to overfitting. Based on this, we enable the estimation of noise
in the training data with Otsu’s algorithm. With this estimation, we train the
network with a modified weighted cross-entropy loss function. We show in our
experiments that Blind Knowledge Distillation detects overfitting effectively
during training and improves the detection of clean and noisy labels on the
recently published CIFAR-N dataset. Code is available at
GitHub111https://github.com/TimoK93/blind_knowledge_distillation.

##  1 Introduction

$0$$0.25$$0.5$$0.75$$1$$0$$0.1$$0.2$$0.3$$0.4$$\mu_{1}$$s$$\mu_{2}$$P(y=\overline{y}|x)$$n$Noisy
LabelsClean LabelsOtsu Distributions Figure 1: Distribution of ground truth
label related probabilities $P_{A}(y=\overline{y}|x)$ at beginning of
overfitting (tipping point) and the resulting gaussian distributions after
Otsu’s algorithm for the dataset Worst. Red bars show the normalized
distribution of noisy labels and green bars of clean labels, respectively.
Note that the gaussian distributions (blue) are scaled for visualization
purposes. Our presented Blind Knowledge Distillation enables an adaptive noise
estimation via the thresholds $\mu_{1}$, $s$, and $\mu_{2}$ and a robust
learning with noisy labels.

Learning with noisy labels is a challenging task in image classification. It
is well known that label noise leads to heavy performance drops with standard
classification methods Song et al. (2022). The goal of learning with noisy
labels is therefore to train a classification model with labelled training
images and achieve high classification performance on unseen test images, even
if the labels for training are noisy and corrupted. Labels are noisy because
humans are naturally unable to classify images perfectly due to ambiguous
images, individual human bias, pressure of time, or various other reasons.
Many modern methods Liu et al. (2022b); Rawat and Wang (2017) are trained on
large and potentially noisy datasets and thus it is an interest of the
community to make classification robust against noisy labels.

To evaluate the robustness of methods for learning with noisy labels, clean
image datasets like CIFAR Krizhevsky and Hinton (2009), Clothing1M Xiao et al.
(2015), or WebVision Li et al. (2017) are synthetically corrupted by randomly
flipping label annotations either symmetrically without constraints or
asymmetrically with predefined rules to mimic realistic label noise. However,
Wei et al. Wei et al. (2022) shows that synthetic label noise has different
behaviour compared to real-world label noise and is thus not an ideal choice
to evaluate robust learning. To close this gap, Wei et al. have made great
efforts and presented CIFAR-N with multiple newly annotated ground truth
labels for CIFAR with human-induced label noise. With these new annotations,
robust learning can be evaluated more realistically.

In this paper, we introduce a novel method to detect the beginning of
overfitting on sample details during training, that is usually roughly
estimated as in Li et al. (2020), and present a simple but effective method to
detect most likely corrupted labels. Our method is inspired by Knowledge
Distillation Hinton et al. (2015) for neural networks which extracts
‘knowledge’ from a teacher network to train a student network. Differently
than usual, our student network is just trained with a subset of the teachers
‘knowledge’. Specifically, it does not ‘see’ the ‘knowledge’ about the given
and potentially corrupted ground truth labels by utilizing the teachers ground
truth complementary logits. Therefore we call it Blind Knowledge Distillation.
Based on the detected noisy labels, we propose a simple but effective loss-
correction method to train the teacher model robustly with label noise. We
perform extensive experiments on CIFAR-10N and the results show that Blind
Knowledge Distillation

  * •

successfully estimates the tipping point from fitting to general patterns to
(over)fitting to sample details,

  * •

is an effective method to estimate the likelihood of labels being noisy,

  * •

and improves the classification accuracy while training with high noise
levels.

##  2 Related Work

Methods in the field of robust learning to tackle noisy labels can be divided
into label correction, loss correction, and refined strategies Song et al.
(2022); Wang et al. (2019). In this section, we contextualize the latest
methods of the CIFAR-N leaderboard based on the aforementioned categories.

Label correction is an approach in which the given ground truth labels are
dynamically changed during optimization to obtain labels of higher quality.
SOP Liu et al. (2022a) performs label correction by optimizing the ground
truth labels with Stochastic Gradient Descent (SGD) w.r.t. the classification
loss. It alternates between the update of model weights and the update of
additional soft-label weights.

Another approach is loss correction which is usually applied by weighting the
loss term or adding a new loss for each sample in the training dataset. The
methods CORES Cheng et al. (2021) and ELR Liu et al. (2020) add a
regularization term to the standard cross-entropy (CE) loss to penalize likely
corrupted labels. PeerLoss Liu and Guo (2020) introduces and minimizes peer
loss functions between randomly selected samples. CAL Zhu et al. (2021)
extends this approach and estimates the covariances between noise rates and
their bayes optimal label.

The last category tackles noisy labels by using refined strategies. CoTeaching
Han et al. (2018) trains a neural network with samples with high confident
predictions of a second network, and vice versa. DivideMix Li et al. (2020)
and PES Bai et al. (2021) split the dataset in clean and corrupted subsets and
apply semi-supervised learning methods. In detail, DivideMix trains two
independent neural networks and splits the set of one network based on the
predictions of the other network to avoid confirmation bias. In contrast to
this, PES applies early stopping of the optimization to every network layer
independently, instead of applying it to the whole network simultaneously, as
usual.

Our method combines a refined strategy to detect most likely corrupted labels
in the first stage and performs loss-correction in the second stage while
incorporating the estimation of likely corrupted labels. While other methods
manually define warm-up epochs, we adapt to the dataset and estimate the
optimal stopping point for the standard CE training. Instead of applying
extensive semi-supervised augmentation methods, we apply a simple sample
dependent loss correction.

##  3 Preliminaries

Given a set of annotated image samples $X$ and a set of classes $C$, the task
of image classification is to assign every sample $x$ from
$(x,\overline{y})\in X$ to the correct class label $y=c\in C$ without prior
knowledge of the correct class label $y$ and a potentially noisy annotation
$\overline{y}$. Modern methods use neural networks $f(\Phi,x)$ to estimate the
probability distribution $P(y=c|x)$ for every class $c$ Liu et al. (2022b); Li
et al. (2020); Liu et al. (2022a); Cheng et al. (2021); He et al. (2016a), in
which $\Phi$ denotes a set of trainable network parameters. More specifically,
neural networks predict a logit vector $\vec{l}\in\mathbb{R}^{|C|}$ with a
logit $l_{c}$ for every class and transform it into probabilities with the
softmax function

$P(y=c|x)=\frac{e^{l_{c}}}{\sum_{i\in C}e^{l_{i}}}.$ (1)

Finally, the class $c$ with the highest probability $P(y=c|x)$ is assumed to
be the correct label $y$.

The task is to define the network architecture of $f$ and the training
strategy to optimize $\Phi$, so that $f(\Phi,x)$ predicts a satisfying
distribution $P(y=c|x)$ in which the correct class has the highest
probability. Most methods optimize $\Phi$ with large manually annotated image
datasets and minimize the categorical cross entropy (CE) loss objective

$L_{\text{CE}}=\frac{1}{|X|}\sum_{\begin{subarray}{c}(x,\overline{y})\in
X\end{subarray}}-\log\big{(}P(y=\overline{y}|x)\big{)}$ (2)

or one of its derivatives.

Extending the task of image classification, the challenging task of learning
with noisy labels addresses the problem that the given ground truth labels
$\overline{y}$ could be noisy and not the true labels $\overline{y}\neq y$.
False ground truth labels dramatically impede the optimization of $\Phi$.
Thus, the goal is to train classifiers with an accuracy that is comparable to
classifiers that would be optimized with clean labels $\overline{y}=y$. A
second goal is to identify noisy labels $\overline{y}\neq y$ in the dataset.

The approach proposed in this paper addresses both tasks. Note that the method
is iteratively trained with random sampled batches $X^{\prime}\subset X$. We
keep the notation of $X$ in the next sections for simplicity, e.g. in Eq. (3).

##  4 Method

Figure 2: Our proposed Blind Knowledge Distillation framework. The teacher and
student network share the same topology but have different weights. While the
teacher network is trained by optimizing the standard CE loss, the student
network is trained with the ground truth complementary logits of the teacher
network and the mean squared error loss. The teacher network predicts the
class membership probabilities $P_{T}(y=c|x)$ and the student network predicts
the probabilities $P_{S}(y=c|x)$.

To enlarge the robustness of neural networks against label noise and to detect
noisy labels, we present a novel training strategy to estimate the likelihood
of every label being noisy and apply a weighted loss based on this.

First, we adapt the student-teacher architecture Gou et al. (2021) and
introduce Blind Knowledge Distillation to extract generalized patterns from
the data. Then, we present a method to detect the beginning of overfitting
with the student network and enable the detection of noisy labels by
estimating four confidence levels of being noisy. Finally, we optimize $\Phi$
with a robust training strategy to train the final classifier. All three steps
are described in the following sections.

###  4.1 Blind Knowledge Distillation

Neural networks with large number of parameters $\Phi$ can memorize the
training examples, so that $P(y=\overline{y}|x)\approx 1$ for every sample $x$
in the training set $X$. Also, neural networks adapt simple patterns during
the early optimization epochs and overfit to specific image details in the
latter ones. As shown in Liu et al. (2020), valid patterns are learned by
maximizing the logit $l_{y=\overline{y}}$ for clean samples in the first
training stages. Subsequently to this early generalization, maximizing the
logits $l_{y\neq\overline{y}}$ of corrupted samples degrades the
classification accuracy. More important for this method is the phenomenon that
the ground truth complementary logits $l_{c\neq\overline{y}}$ are also
minimized in the latter stages.

To avoid the maximization of $l_{y\neq\overline{y}}$ for noisy labels, we
create a new student-teacher architecture, in which the student only learns
generalized patterns. In the student-teacher architecture, a student model is
trained with the output of a teacher model. This method is called Knowledge
Distillation Hinton et al. (2015) and transfers the patterns that are encoded
in the teacher model to the student. Model bias or wrong patterns from the
teacher can also be transferred. To avoid this undesired transfer, we
introduce the ground truth annotation complementary ‘knowledge’ by removing
all information that is immediately connected to a potentially corrupted
ground truth label $\overline{y}$.

Unlike the usual knowledge distillation architecture, our models share the
same topology $f(\Phi,x)$ but have different weights $\Phi_{T}$ (teacher) and
$\Phi_{S}$ (student). The teacher model $f(\Phi_{T},x)$ is trained with the
standard CE loss (Eq. (2)). The student model is trained with the unlabelled
ground truth complementary logits $l^{T}_{c\neq\overline{y}}$ derived from the
teacher model and an extended but simple mean squared error loss

$L_{\text{Stud}}=\frac{1}{|X|}\sum_{\begin{subarray}{c}(x,\overline{y})\in
X\end{subarray}}\frac{1}{|C|-1}\sum_{\begin{subarray}{c}c\in C\\\
c\neq\overline{y}\end{subarray}}(l^{T}_{c}-l^{S}_{c})^{2}$ (3)

in which $l^{S}_{c}$ denotes the logits of the student model. This loss
function transfers the generalized patterns by imitating the output of the
teacher model but without taking the potentially corrupted ground truth label
$\overline{y}$ into account. The training architecture is visualized in Fig.
2.

Since high valued complementary logits $l_{c\neq\overline{y}}$ are minimized
in the latter optimization stages after general features are learned, also the
resulting probabilities after softmax (Eq. (1)) converge to a uniform
distribution. Thus, we can identify approximately the training epoch, in which
the neural network starts overfitting to specific sample details by monitoring
the mean maximal probability

$\hat{p}_{\text{max}}=\frac{1}{|X|}\sum_{x\in X}\max_{c\in
C}\big{(}P(y=c|x)\big{)}$ (4)

of the students network. During training in epoch $i$, the fitting-epoch in
which $\hat{p}^{i}_{\text{max}}$ is maximal can be certainly identified online
with a delay of $k$ epochs by checking if $\hat{p}^{i-k}_{\text{max}}$ is the
maximum of the last $2k+1$ epochs.

Furthermore, it shows that the student model has the ability to classify
images comparable to the teacher model before the detail fitting starts. Thus,
we modify our final classification probability by combining the teacher’s
prediction $P_{T}(y=c|x)$ and the student’s prediction $P_{S}(y=c|x)$ to the
agreement probability

$P_{A}(y=c|x)=\frac{P_{T}(y=c|x)\cdot P_{S}(y=c|x)}{\sum_{i\in
C}P_{T}(y=i|x)\cdot P_{S}(y=i|x)}.$ (5)

$P_{A}$ enables the noise estimation in the dataset described in the
following. The classification accuracy $P_{T}$, $P_{S}$, and $P_{A}$ are
elaborated more in detail in the experiments (see Sec. 5.4).

###  4.2 Adaptive Noise Estimation

The knowledge about the presence of noise can be used to apply loss
correction. Unfortunately, this knowledge is not given, so we estimate the
probability of a data sample to be noisy. We split the dataset into four
subsets based on Otsu’s algorithm Otsu (1979), in which the membership to a
subset indicates the likelihood of being noisy. Given a set of data samples
$(x,\overline{y})\in X$ with their corresponding agreement probability
$P_{A}(y=\overline{y}|x)$, the first step is to find a threshold $s$ that
splits $X$ into two distributions $X_{1}=\\{x\in X|P_{A}(y=\overline{y}|x)\leq
s\\}$ and $X_{\text{2}}=\\{x\in X|P_{A}(y=\overline{y}|x)>s\\}$, in which
$X_{1}$ contains images with likely noisy and $X_{\text{2}}$ images with
likely clean labels. To find $s$, we assume that $X_{1}$ and $X_{2}$ can be
approximated by two gaussian distributions $(\mu_{1},\sigma_{1})$ and
$(\mu_{2},\sigma_{2})$. The optimal threshold $s$ maximizes the objective

$Q(s)=\frac{n_{1}(s)\big{(}\mu_{1}(s)-\mu\big{)}^{2}+n_{2}(s)\big{(}\mu_{2}(s)-\mu\big{)}^{2}}{n_{1}(s)\sigma_{1}(s)^{2}+n_{2}(s)\sigma_{2}(s)^{2}},$
(6)

where $n_{1}(s)$ and $n_{2}(s)$ denote the cardinality of $X_{1}$ and $X_{2}$
depending on $s$, and $\mu$ is the mean probability $P_{A}(y=\overline{y}|x)$
of all samples in $X$. The optimal $s$ minimizes the inter-class variance and
can be found by calculating $Q(s)$ for all $s$ with a reasonable step size
$\Delta s=0.001$.

Using Otsu’s algorithm, we preserve a threshold $s$ to split the data into
noisy and clean samples, and furthermore thresholds $\mu_{1}$ and $\mu_{2}$ to
subdivide the subsets into more fine-grained subsets. A finer distinction
w.r.t. the likelihood of being noisy allows a more precise weighting of the
samples in the following steps. Depending on the requirements of the
application, the task of label noise detection can be solved by classifying a
sample $x$ by comparing $P_{A}(y=\overline{y}|x)$ with one of the thresholds
$s$, $\mu_{1}$, and $\mu_{2}$. While using $\mu_{1}$ is more liberal to
classifying noisy labels into the clean dataset than $s$, $\mu_{2}$ is more
conservative. A visualization of a distribution of $P_{A}(y=\overline{y}|x)$
and the estimated noise is shown in Fig. 1.

###  4.3 Robust Optimization

After splitting up the dataset into potentially clean and corrupted data, we
use simple robust training techniques to train the final classification model.
Based on the idea of label smoothing Szegedy et al. (2016), we extend the CE
loss (Eq. (2)) and combine the ground truth label with the student’s
prediction and a sample dependent $\alpha_{x}$:

$L_{\text{Robust}}=-\frac{1}{|X|}\sum_{(x,\overline{y})\in X}\sum_{c\in
C}\mathcal{S}(\beta^{c}_{x})\log\big{(}P_{T}(y=c|x)\big{)}\\\ \text{with}\
\beta^{c}_{x}=(1-\alpha_{x})\mathbbm{1}[c=\overline{y}]+\alpha_{x}P_{S}(y=c|x)$
(7)

and Sharpening $\mathcal{S}$ that is explained later.

While the teacher network is trained by $L_{\text{Robust}}$, the student
network is still trained with the student loss $L_{\text{Stud}}$ (Eq. (3)). As
larger as the instance dependent $\alpha_{x}$ gets, the less the ground truth
of a sample $x$ is trusted. We adapt $\alpha_{x}$ for every sample
individually, depending on the cluster membership after Otsu. We define four
fixed alpha values with $\alpha^{1}<\alpha^{2}<\alpha^{3}<\alpha^{4}$ where
$\alpha^{1}$ gets assigned to samples with
$P_{A}(y=\overline{y}|x)\geq\mu_{2}$, $\alpha^{2}$ to samples with
$\mu_{2}>P_{A}(y=\overline{y}|x)\geq s$, $\alpha^{3}$ to samples with
$s>P_{A}(y=\overline{y}|x)\geq\mu_{1}$, and $\alpha^{4}$ otherwise.

Since a larger $\alpha_{x}$ enlarges the entropy in the objective, we use a
modified Sharpening method

$\mathcal{S}(\beta^{c}_{x})=\frac{(\beta^{c}_{x})^{1+\alpha_{x}}}{\sum_{i\in
C}{(\beta^{i}_{x})^{1+\alpha_{x}}}}$ (8)

as used by Li et al. (2020) to minimize the entropy. The Sharpening function
is applied stronger for insecure samples by reusing the above mentioned alpha.

##  5 Experiments

Acc [%] | Aggre | Rand1 | Rand2 | Rand3 | Worst
---|---|---|---|---|---
SOP | 95.61 | 95.28 | 95.31 | 95.39 | 93.24
CORES | 95.25 | 94.45 | 94.88 | 94.47 | 91.66
DivideMix | 95.01 | 95.16 | 95.23 | 95.21 | 92.56
ELR+ | 94.83 | 94.43 | 94.20 | 94.34 | 91.09
PES | 94.66 | 95.06 | 95.19 | 95.22 | 92.68
ELR | 92.38 | 91.46 | 91.61 | 91.41 | 83.58
CAL | 91.97 | 90.93 | 90.75 | 90.74 | 85.36
CE | 87.77 | 85.02 | 86.14 | 85.16 | 77.69
Ours | 93.68 | 92.50 | 92.63 | 92.54 | 86.64

Table 1: Classification accuracy of our method compared to standard CE-loss framework and state-of-the-art methods SOP Liu et al. (2022a), CORES Cheng et al. (2021), DivideMix Li et al. (2020), PES Bai et al. (2021), ELR Liu et al. (2020), and CAL Zhu et al. (2021) . $[\%]$ | Aggre | Rand1 | Rand2 | Rand3 | Worst
---|---|---|---|---|---
$F_{1}$ | Pr | Re | $F_{1}$ | Pr | Re | $F_{1}$ | Pr | Re | $F_{1}$ | Pr | Re | $F_{1}$ | Pr | Re
$P_{T}$ |  $\mu_{1}$ | 69.2 | 58.0 | 85.9 | 82.9 | 82.7 | 83.0 | 83.4 | 83.6 | 83.3 | 83.1 | 82.9 | 83.2 | 75.0 | 95.9 | 61.6
$s$ | 46.6 | 30.7 | 97.1 | 67.2 | 51.3 | 97.5 | 67.9 | 52.1 | 97.6 | 67.3 | 51.3 | 97.6 | 85.9 | 79.6 | 93.3
$\mu_{2}$ | 28.8 | 16.8 | 99.7 | 44.8 | 28.8 | 99.8 | 45.9 | 29.8 | 99.8 | 45.1 | 29.1 | 99.8 | 70.8 | 55.0 | 99.6
$P_{S}$ |  $\mu_{1}$ | 50.9 | 35.7 | 88.4 | 74.5 | 65.7 | 86.0 | 75.8 | 67.9 | 85.9 | 74.8 | 66.1 | 86.0 | 75.8 | 92.3 | 64.4
$s$ | 31.9 | 19.0 | 98.3 | 55.0 | 38.2 | 98.2 | 56.6 | 39.7 | 98.1 | 55.4 | 38.6 | 98.5 | 81.9 | 72.2 | 94.4
$\mu_{2}$ | 22.3 | 12.5 | 99.8 | 39.8 | 24.8 | 99.9 | 41.3 | 26.0 | 99.9 | 40.4 | 25.3 | 99.9 | 68.9 | 52.6 | 99.7
$P_{A}$ |  $\mu_{1}$ | 71.2 | 63.7 | 80.5 | 82.5 | 84.1 | 81.0 | 83.1 | 85.4 | 80.9 | 82.7 | 84.9 | 80.7 | 77.8 | 95.8 | 65.5
$s$ | 55.7 | 39.6 | 94.1 | 73.7 | 60.1 | 95.2 | 75.2 | 62.1 | 95.3 | 74.6 | 61.3 | 95.1 | 87.3 | 83.4 | 91.5
$\mu_{2}$ | 37.0 | 22.8 | 98.9 | 54.2 | 37.3 | 99.4 | 55.4 | 38.4 | 99.2 | 54.9 | 37.9 | 99.4 | 75.7 | 61.3 | 99.0
Table 2: $F_{1}$-Score, Precision, and Recall on the label noise detection
task with different probability sets ($P_{T}$, $P_{S}$, $P_{A}$) and different
thresholds provided by Otsu’s method ($s$, $\mu_{1}$, $\mu_{2}$) to split the
dataset into clean and corrupted subsets. Best metrics are presented in bold.

$0$$25$$50$$75$$0.2$$0.3$$0.4$$0.5$$0.6$Epoch$\hat{p}_{\text{max}}$Average
Maximal Prediction (Student)WorstRand1Aggre (a) The average maximal
probability of the students network prediction. Vertical lines denote the
maxima during the training process which can be interpreted as tipping point
at which the model start fitting to individual sample details.

$0$$25$$50$$75$$0.2$$0.4$$0.6$$0.8$EpochAccClassification Accuracy
(Teacher)WorstRand1Aggre (b) Test accuracy of the teacher network during
training. Vertical lines denote the tipping points from Fig. 3. With
increasing noise rates, the detected tipping points fit approximately to the
maxima of the classification test accuracy.

$0$$25$$50$$75$$0$$0.2$$0.4$$0.6$$0.8$Epoch$F_{1}$Detection $F_{1}$-Score
(Agreement)WorstRand1Aggre (c) Noise detection $F_{1}$-score on the training
set with student maxima from Fig. 3. Dashed lines indicate split by Otsu
threshold $s$, while solid lines indicate the split by $\mu_{1}$. The student
maxima fit approximately to the maxima of the detection performance for
$\mu_{1}$ for higher noise rates.

Figure 3: Relation between the student’s prediction behavior and the
classification and noise detection performance. The maximum of the student
network’s average maximum probability indicates the start of fitting to sample
details and thus can be used for early stopping to avoid overfitting. The
models for this figure are trained without robust training for 75 epochs to
show the standard training behavior. $0$$25$$50$$75$$0.6$$0.8$EpochAcc$P_{T}$
(Teacher)$P_{S}$ (Student)$P_{A}$ (Agreement) Figure 4: Test accuracy w.r.t.
training epoch based on the predicted probabilities of the teacher, student,
and the proposed combined agreement. The agreement probability combines the
strengths and outperforms the teacher and student probability. Note that we
train the framework for 75 epochs and without robust optimization (Sec. 4.3).

We perform several experiments to evaluate our proposed method. The
experimental setup and the used metrics are explained first. Then we present
evaluation metrics on the recently released dataset CIFAR-10N Wei et al.
(2022) and show details and observations of our core method Blind Knowledge
Distillation.

###  5.1 Experimental Setup

We evaluate our method on the noise levels provided in the CIFAR-10N dataset.
To be comparable to other methods, we utilize the same model setup as used in
Li et al. (2020). We use a 18-layer PreAct ResNet He et al. (2016b) and
Stochastic Gradient Descent with momentum of $0.9$ and weight decay of
$0.0005$ as optimizer. The networks are trained for $300$ epochs beginning
with a learning rate of $0.02$ and reduce it to $0.002$ after $150$ epochs. We
train the network with randomly sampled batches of $128$ image samples. In the
first stage, the teacher network optimizes the standard CE loss (Eq. (2))
until the detection of the tipping point induces the start of the second
stage, in which the teacher network optimizes the modified loss (Eq. (7)). The
hyperparameters introduced by our method are set to $\alpha^{1}=0.3$,
$\alpha^{2}=0.45$, $\alpha^{3}=0.55$, $\alpha^{4}=0.7$. For the noisy
detection task, we use the probabilities $P_{A}$ and the threshold $\mu_{1}$
after Otsu. We repeated the experiments at least five times with random seeds
and report the averaged metrics.

The method proposed in this paper is evaluated on CIFAR-10N Wei et al. (2022).
CIFAR-10N manually re-labelled the CIFAR-10 Krizhevsky and Hinton (2009) by
multiple humans to investigate the impact of realistic label noise compared to
synthetically induced ones. The dataset contains 50K training images and 10K
test images with a size of $32\times 32$. For the training set, there are five
label sets with realistic human label noise with a ratio of approx. $9\%$,
$17\%$, $18\%$, $18\%$, and $40\%$ label noise. In the same order of the noise
ratios, we denote them as Aggre, Rand1, Rand2, Rand3, and Worst in our
experiments.

The tasks for the dataset are twofold: First, the classifier should be trained
robust to achieve a high test accuracy even with high noisy rates and second,
noisy labels in the training data should be detected and marked as noisy. The
metrics to evaluate the tasks are given in the next section.

###  5.2 Metrics

We evaluate the performance of image classification with the commonly used
Accuracy (Acc) metric. It measures the classification accuracy on the test
dataset $X_{\text{Test}}$ using the ratio of correct classified test samples
compared to the dataset size:

$\text{Acc}=\frac{\sum_{(x,\overline{y})\in
X_{\text{Test}}}\mathbbm{1}\big{[}\operatorname*{arg\,max}_{c\in
C}\big{(}P(y=c|x)\big{)}=\overline{y}\big{]}}{|X_{\text{Test}}|}$ (9)

The task of noisy label detection is evaluated with the well-known
$F_{1}$-score, Precision (Pr), and Recall (Re) metrics, in which Precision
decreases if clean labels are classified as noisy and Recall decreases if
noisy labels are classified as clean. The $F_{1}$-score harmonizes both
aspects. With the subsets of true ($X_{\text{Noise}}\subset X$) and predicted
($X^{\prime}_{\text{Noise}}\subset X$) noisy labels from the training set $X$,
the metrics are defined as:

$\text{Pr}=\frac{\sum_{x\in X^{\prime}_{\text{Noise}}}\mathbbm{1}\big{[}x\in
X_{\text{Noise}}]}{|X^{\prime}_{\text{Noise}}|},$ (10)
$\text{Re}=\frac{\sum_{x\in X_{\text{Noise}}}\mathbbm{1}\big{[}x\in
X^{\prime}_{\text{Noise}}]}{|X_{\text{Noise}}|},$ (11) $\text{and}\quad
F_{1}=\frac{2}{\text{Pr}^{-1}+\text{Re}^{-1}}.$ (12)

###  5.3 CIFAR-10N

This section elaborates the results for the tasks of robust training and noise
detection.

#### Robust Training

We compare our results to the latest six state-of-the-art methods and the
standard CE baseline on the CIFAR-10N Leaderboard in Tab. 1. Our method
achieves the performance to be listed on the new sixth position of the
leaderboard outperforming CAL and standard ELR. We want to mention that ELR+
and DivideMix apply multiple models and high performance semi-supervised
strategies such as MixMatch Berthelot et al. (2019).

#### Noise Detection

The detection performance is shown in Tab. 2. We present $F_{1}$, precision,
and recall for all five noise levels in CIFAR-10N (10 classes). The split to
classify clean and corrupted labels is performed based on one of the
probability sets $P_{T}$, $P_{S}$, and $P_{A}$ and the three threshold $s$,
$\mu_{1}$, $\mu_{2}$ provided by Otsu’s method. Intuitively, the precision is
higher for the lower threshold $\mu_{1}$ and recall for the higher threshold
$\mu_{2}$, respectively. The experiments show that the harmonized metric
$F_{1}$ performs best for $\mu_{1}$. Thus, $\mu_{1}$ is used to solve the task
of Noise detection. The combined probability $P_{A}$ performs better or on par
w.r.t. the $F_{1}$-score, confirming the improved classification accuracy that
is also visible in Fig. 3. An exemplary distribution of ground truth label
probabilities for clean and noisy labels with the subsequent split based on
Otsu is shown for the Worst dataset in Fig. 1.

###  5.4 Blind Knowledge Distillation

The core contribution of our method is Blind Knowledge Distillation. This
section analyzes its ability to detect (over)fitting on sample details and how
it can be used to improve the classification accuracy. For the experiments in
this section, we trained our student-teacher framework without detection of
the tipping point and robust optimization (Sec. 4.3) after generalization for
75 epochs to show the teacher’s and student’s learning behavior.

Fig. 3 shows the average maximal probability prediction of the student network
over training time. The probability strongly increases in the first training
epochs and degrades after a tipping point. It is notable that high noise rates
decrease the absolute mean probability in general (see Worst). Our explanation
of this behavior is that classifiers are not able to clearly predict a class
based on simple and generalized but ambiguous image patterns. Thus, the
classifier produces multiple predictions $P(y=c|x)\gg 0$ during the first
generalization stage. An example pattern could be the coarse shape which is
often ambiguous, e.g. for classes dog and cat. In the second training stage
after the tipping point, the teacher network adapts detailed sample patterns
to maximize $P(y=\overline{y}|x)$ which also minimizes
$P(y\neq\overline{y}|x)$

The test classification accuracy of the teacher network is shown in Fig. 3.
While early stopping of standard optimization is not important for clean
datasets or low noise levels (Aggre), fitting on sample details leads to
overfitting and decreases the classification accuracy on high noise rates
(Worst). Therefore, the choice of an early stopping epoch is highly important.
It shows that the tipping point from Fig. 4 is a good indicator to detect
overfitting. It proposes an accurate estimation to stop optimizing on high
noise levels without stopping too early on low noise levels.

Using the tipping point in Fig. 3 to split the data into potentially corrupted
and clean subsets is intuitively, due to the beginning overfitting and
decreasing classification accuracy. Fig. 3 shows the detection ability with
the $F_{1}$-score if splitting the dataset based on Otsu’s algorithm and
$P_{A}$ at every epoch. Similar to Fig. 3, the tipping point gives a guess for
a suitable epoch to split the dataset. The estimation for high noise rates is
sufficient, especially due to the decreasing $F_{1}$-score after the tipping
point. Since low noise rates does not seem to affect the $F_{1}$ negatively in
latter training stages, the tipping point estimate leads to a slightly too
early splitting epoch.

An interesting insight about the teacher and student classification accuracy
$P_{T}$ and $P_{S}$ is shown in Fig. 4 on a high noise level. While the
teachers accuracy decreases during overfitting, the students accuracy
persists. We claim that using the complementary student loss from Eq. (3)
prevents the student from fitting to misleading image details by removing the
ground truth related logits $l_{c=\overline{y}}$. Also interesting is that the
combined probability $P_{A}$ suits as the overall best probability for
classification. While the combined probability is quite similar to the
teachers probability $P_{A}\approx P_{T}$ during the early stage, it converges
to the students accuracy $P_{A}\approx P_{S}$ in the latter ones. Near the
tipping point, it outperforms both.

Overall, Blind Knowledge Distillation is a better choice to automatically
detect overfitting rather than to stop training after predefined and fixed
periods (e.g. in Li et al. (2020)). Combining $P_{T}$ and $P_{S}$ to $P_{A}$
can be used to improve the overall classification accuracy.

##  6 Conclusion

This paper introduces Blind Knowledge Distillation that is able to transfer
simple and general image patterns that are not based on individual image
details. We show that our framework is able to identify the tipping point from
fitting to simple but general image patterns to fitting to image details and
use it for early stopping in standard classification frameworks and
furthermore to estimate the likelihood of samples in the training data of
being clean or corrupted.

Our method performs on par with state-of-the-art methods that are not extended
with high performance semi-supervised training strategies. Compared to them,
we do not rely on manually predefined warm-up phases and adapt it online
during training. However, the intention of this paper is to provide new
insights about general learning behavior rather than to tune our method with
known strategies. We hope that Blind Knowledge Distillation helps researchers
to improve the handling of under- and overfitting.

## Acknowledgments

This work was supported by the Federal Ministry of Education and Research
(BMBF), Germany under the project LeibnizKILabor (grant no. 01DD20003) and the
Deutsche Forschungsgemeinschaft (DFG) under Germany’s Excellence Strategy
within the Cluster of Excellence PhoenixD (EXC 2122).

## References

  * Bai et al. [2021] Yingbin Bai, Erkun Yang, Bo Han, Yanhua Yang, Jiatong Li, Yinian Mao, Gang Niu, and Tongliang Liu.  Understanding and improving early stopping for learning with noisy labels.  In M. Ranzato, A. Beygelzimer, Y. Dauphin, P.S. Liang, and J. Wortman Vaughan, editors, Advances in Neural Information Processing Systems, volume 34. Curran Associates, Inc., 2021.
  * Berthelot et al. [2019] David Berthelot, Nicholas Carlini, Ian Goodfellow, Nicolas Papernot, Avital Oliver, and Colin A Raffel.  Mixmatch: A holistic approach to semi-supervised learning.  In H. Wallach, H. Larochelle, A. Beygelzimer, F. d'Alché-Buc, E. Fox, and R. Garnett, editors, Advances in Neural Information Processing Systems, volume 32. Curran Associates, Inc., 2019.
  * Cheng et al. [2021] Hao Cheng, Zhaowei Zhu, Xingyu Li, Yifei Gong, Xing Sun, and Yang Liu.  Learning with instance-dependent label noise: A sample sieve approach.  In International Conference on Learning Representations, 2021.
  * Gou et al. [2021] Jianping Gou, Baosheng Yu, Stephen J Maybank, and Dacheng Tao.  Knowledge distillation: A survey.  International Journal of Computer Vision, 129(6), 2021.
  * Han et al. [2018] Bo Han, Quanming Yao, Xingrui Yu, Gang Niu, Miao Xu, Weihua Hu, Ivor Tsang, and Masashi Sugiyama.  Co-teaching: Robust training of deep neural networks with extremely noisy labels.  In S. Bengio, H. Wallach, H. Larochelle, K. Grauman, N. Cesa-Bianchi, and R. Garnett, editors, Advances in Neural Information Processing Systems, volume 31. Curran Associates, Inc., 2018.
  * He et al. [2016a] Kaiming He, Xiangyu Zhang, Shaoqing Ren, and Jian Sun.  Deep residual learning for image recognition.  In Proceedings of the IEEE conference on computer vision and pattern recognition, 2016.
  * He et al. [2016b] Kaiming He, Xiangyu Zhang, Shaoqing Ren, and Jian Sun.  Identity mappings in deep residual networks.  In European conference on computer vision. Springer, 2016.
  * Hinton et al. [2015] Geoffrey Hinton, Oriol Vinyals, and Jeff Dean.  Distilling the knowledge in a neural network.  arXiv preprint arXiv:1503.02531, 2015.
  * Krizhevsky and Hinton [2009] Alex Krizhevsky and Geoffrey Hinton.  Learning multiple layers of features from tiny images.  2009\.
  * Li et al. [2017] Wen Li, Limin Wang, Wei Li, Eirikur Agustsson, and Luc Van Gool.  Webvision database: Visual learning and understanding from web data.  arXiv preprint arXiv:1708.02862, 2017.
  * Li et al. [2020] Junnan Li, Richard Socher, and Steven C.H. Hoi.  Dividemix: Learning with noisy labels as semi-supervised learning.  In International Conference on Learning Representations, 2020.
  * Liu and Guo [2020] Yang Liu and Hongyi Guo.  Peer loss functions: Learning from noisy labels without knowing noise rates.  In International Conference on Machine Learning. PMLR, 2020.
  * Liu et al. [2020] Sheng Liu, Jonathan Niles-Weed, Narges Razavian, and Carlos Fernandez-Granda.  Early-learning regularization prevents memorization of noisy labels.  In H. Larochelle, M. Ranzato, R. Hadsell, M.F. Balcan, and H. Lin, editors, Advances in Neural Information Processing Systems, volume 33. Curran Associates, Inc., 2020.
  * Liu et al. [2022a] Sheng Liu, Zhihui Zhu, Qing Qu, and Chong You.  Robust training under label noise by over-parameterization.  arXiv preprint arXiv:2202.14026, 2022.
  * Liu et al. [2022b] Ze Liu, Han Hu, Yutong Lin, Zhuliang Yao, Zhenda Xie, Yixuan Wei, Jia Ning, Yue Cao, Zheng Zhang, Li Dong, Furu Wei, and Baining Guo.  Swin transformer v2: Scaling up capacity and resolution.  In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, 2022.
  * Otsu [1979] Nobuyuki Otsu.  A threshold selection method from gray-level histograms.  IEEE transactions on systems, man, and cybernetics, 9(1), 1979.
  * Rawat and Wang [2017] Waseem Rawat and Zenghui Wang.  Deep convolutional neural networks for image classification: A comprehensive review.  Neural Computation, 29(9), 2017.
  * Song et al. [2022] Hwanjun Song, Minseok Kim, Dongmin Park, Yooju Shin, and Jae-Gil Lee.  Learning from noisy labels with deep neural networks: A survey.  IEEE Transactions on Neural Networks and Learning Systems, 2022\.
  * Szegedy et al. [2016] Christian Szegedy, Vincent Vanhoucke, Sergey Ioffe, Jon Shlens, and Zbigniew Wojna.  Rethinking the inception architecture for computer vision.  In 2016 IEEE Conference on Computer Vision and Pattern Recognition, 2016.
  * Wang et al. [2019] Yisen Wang, Xingjun Ma, Zaiyi Chen, Yuan Luo, Jinfeng Yi, and James Bailey.  Symmetric cross entropy for robust learning with noisy labels.  In Proceedings of the IEEE/CVF International Conference on Computer Vision (ICCV), October 2019.
  * Wei et al. [2022] Jiaheng Wei, Zhaowei Zhu, Hao Cheng, Tongliang Liu, Gang Niu, and Yang Liu.  Learning with noisy labels revisited: A study using real-world human annotations.  In International Conference on Learning Representations, 2022.
  * Xiao et al. [2015] Tong Xiao, Tian Xia, Yi Yang, Chang Huang, and Xiaogang Wang.  Learning from massive noisy labeled data for image classification.  In 20155 IEEE Conference on Computer Vision and Pattern Recognition, 2015.
  * Zhu et al. [2021] Zhaowei Zhu, Tongliang Liu, and Yang Liu.  A second-order approach to learning with instance-dependent label noise.  In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, 2021.


# Phase tuning of multiple Andreev reflections of Dirac fermions and the
Josephson supercurrent in Al-MoTe2-Al junctions.

Zheyi Zhu1 Stephan Kim1 Shiming Lei2 Leslie M. Schoop2 R. J. Cava2 N. P.
Ong1,§ 1Department of Physics, Princeton University, Princeton, NJ 08544, USA
2Department of Chemistry, Princeton University, Princeton, NJ 08544, USA

###### Abstract

When a normal metal $N$ is sandwiched between two superconductors, the energy
gaps in the latter act as walls that confine electrons in $N$ in a square-well
potential. If the voltage $V$ across $N$ is finite, an electron injected into
the well undergoes multiple Andreev reflections (MAR) until it gains enough
energy to overcome the energy barrier. Because each reflection converts an
electron to a hole (or vice versa), while creating (or destroying) a Cooper
pair, the MAR process shuttles a stream of pairs across the junction. An
interesting question is, given a finite $V$, what percentage of the shuttled
pairs end up as a Josephson supercurrent? This fraction does not seem to have
been measured. Here we show that, in high-transparency junctions based on the
type II Dirac semimetal MoTe2, the MAR leads to a stair-case profile in the
current-voltage ($I$-$V$) response, corresponding to pairs shuttled
incoherently by the $n^{th}$-order process. By varying the phase $\varphi$
across the junction, we demonstrate that a Josephson supercurrent ${\bf
J}_{\rm s}\sim\sin\varphi$ co-exists with the MAR steps, even at large $V$.
The observed linear increase in the amplitude of ${\bf J}_{\rm s}$ with $n$
(for small $n$) implies that ${\bf J}_{\rm s}$ originates from the population
of pairs that are coherently shuttled. We infer that the MAR steps and the
supercurrent are complementary aspects of the Andreev process. The experiment
yields the percentage of shuttled pairs that form the supercurrent. At large
$V$, the coherent fraction is initially linear in $n$. However, as $V\to 0$
($n\gg 1$), almost all the pairs end up as the observed Josephson
supercurrent.

An active area of research is the investigation of proximity-induced pairing
correlations in topological and other unconventional systems. Examples are
graphene Herrero ; DuAndrei , carbon nanotubes Pillet , point-contact or break
junctions Esteve ; Pothier , Josephson $\omega$-junctions Strambini ,
topological Bi nanowires Bouchiat , and systems exhibiting edge currents
Yacoby ; Kim . In these experiments, Andreev reflections and the associated
subgap or bound states play central roles. Signatures specific to topological
junctions arising from Andreev subgap states have been discussed by several
groups SanJose ; vonOppen ; Halperin .

In the MAR process (Fig. 1a), a right-moving electron in $N$ (red circle) is
Andreev reflected at the right $N$-$S$ interface as a left-moving hole
(white), which is, in turn, reflected as an electron at the left interface. If
the voltage $V$ is finite, both excitations gain energy $eV$ with each
traversal ($e>0$ is the elemental charge). Eventually, with $n$ traversals
(the $n^{\rm th}$-order process), the excitation acquires enough energy to
surmount the potential barrier. Two successive Andreev reflections shuttle one
Cooper pair across the junction (green arrows). Previously, calculations have
shown that subgap states in point-contact junctions give rise to a
supercurrent when $V=0$ Furusaki ; Beenakker ; Furusaki99 . We provide
evidence that, even at large $V$, a supercurrent shuttled by the MAR process
persists in large $S$-$N$-$S$ junctions.

Weak subgap, subharmonic features have long been observed in numerous
experiments on (single) $S$-${\cal I}$-$S$ junctions (${\cal I}$ = insulator)
and ascribed to various causes Taylor ; Marcus . Their identification with MAR
was made in Ref. Klapwijk . Subsequently, microscopic calculations of MAR were
compared with experiments on point-contact or break junctions Averin ; Bratus
; Cuevas . Here we extend these pioneering experiments to a new regime,
applying the powerful technique of phase tuning on high-transparency $SNS$
junctions using an asymmetric SQUID (superconducting quantum interference
device) layout Ouboter ; Fulton ; Barone ; Esteve ; Pothier ; Bouchiat ;
vonOppen . Our focus is on the MAR and the Josephson effect when $N$ is the
Dirac-Weyl semimetal MoTe2 (we stay above its critical temperature $T_{\rm c}$
= 100 mK Kim ).

Flux-grown crystals of MoTe2 exhibiting high residual resistivity ratios
($\sim$1 000) and mobilities (100 000 to 150 000 cm2/Vs) were exfoliated in
air into thin flakes of thickness 100 nm, and transfered onto a silicon
substrate capped with a 90 nm-thick SiO2 layer. We deposited Al wires on the
surface of a flake to define four DC SQUIDs (details in Methods). As shown in
Fig. 1b, inset, the (sample) $S$-$N$-$S$ junction (1), with critical current
$I_{\rm c}$, is fabricated in parallel with an auxiliary conventional
$S$-${\cal I}$-$S$ tunnel junction (0) with critical current $I_{0}$ (where
$S$ = Al, $N$ = MoTe2 and ${\cal I}$ = Al2O3). The phases across the
$S$-$N$-$S$ and $S$-${\cal I}$-$S$ junctions are $\delta$ and $\gamma$,
respectively. In the 4 devices (1–4), the ratio $I_{0}/I_{\rm c}$ is $\geq 4$.
The junction width $d$ equals 200, 300, 400 and 500 nm in devices S1, S2, S3
and S4, respectively.

When the applied current $I$ exceeds the SQUID’s critical current $I_{\rm
cS}$, the voltage $V$ rises steeply. We recorded both the $I$-$V$ curve and
the differential resistance $dV/dI$ vs. $I$ in a magnetic field $B$ at
temperatures $T$ from 0.135 to 1.1 K. The curves of $V(B,I)$ and $dV/dI(B,I)$
are reported as color maps in the $B$-$I$ plane (the two experimentally
controlled quantities).

In the regime $V=0$ ($I<I_{\rm cS}$), both phases $\delta$ and $\gamma$ are
static and related by the constraint $\delta-\gamma=\varphi$, where
$\varphi\equiv 2\pi BA/\phi_{0}$ is the flux-induced phase shift and
$\phi_{0}=h/2e$ is the superconducting flux quantum ($A$ is the loop area and
$h$ is Planck’s constant). For $V=0$, we maximize $I_{\rm cS}$ under the
constraint on $\delta-\gamma$ and find that the auxiliary phase $\gamma$ is
pinned close to $\pi/2$ if $I_{0}\gg I_{\rm c}$. The curve of $I_{\rm cS}$ vs.
$\varphi$ then yields the CPR (current-phase-relation) curve of the SNS
junction. From the CPR, we obtain $I_{0}=93\,\mu{\rm A}$ and $I_{\rm
c}=29.6\,\mu$A (see Fig. S1 in Methods).

Our focus is on the finite-$V$ regime ($I>I_{\rm cS}$) Ouboter . At finite
$V$, $\delta(t)$ and $\gamma(t)$ wind rapidly at the same rate
($\dot{\delta}=\dot{\gamma}$) with their difference fixed at $\varphi$. The
finite voltage across the SQUID, given by
$V=\hbar\langle\dot{\gamma}\rangle/2e$ ($\langle\cdots\rangle$ denotes time-
averaged), drives a normal current $I_{\rm N}$ that flows parallel to the
supercurrents $I_{s0}$ and $I_{s1}$ in the auxliary and sample junctions,
respectively ($I=I_{s0}+I_{s1}+I_{\rm N}$).

To gain insight into the MAR, it is helpful to generalize the RSJ model. We
assume that, at finite $V$, the supercurrent in the $S$-$N$-$S$ junction can
be approximated by $I_{s1}\simeq I_{\rm A}(V)\sin\delta$, with a $V$-dependent
amplitude $I_{\rm A}(V)$ (as $V\to 0$, $I_{\rm A}(0)=I_{\rm c}$). The total
$t$-dependent current is then

$I(V,\varphi,t)=\frac{\hbar\dot{\gamma}(t)}{2eR_{\parallel}}+I_{0}\sin\gamma(t)+I_{\rm
A}(V)\sin(\gamma(t)+\varphi),$ (1)

with $R_{\parallel}^{-1}=G=G_{0}+G_{1}$ where $G_{1}$ is the shunt $S$-$N$-$S$
conductance arising from the MAR and $G_{0}$ is the remaining background
conductance.

We consider the two extremal cases when $\varphi$ equals $\varphi_{+}$ (with
$I_{s1}$ parallel to $I_{s2}$) and $\varphi_{-}$ (antiparallel). At
$\varphi_{\pm}$, we have

$I(V,\varphi_{\pm},t)=\frac{\hbar\dot{\gamma}(t)}{2eR_{\parallel}}+[I_{0}\pm
I_{\rm A}(V)]\sin\gamma(t),$ (2)

which reduces to the equation governing the phase dynamics in a single
junction Ouboter ; Fulton ; Barone . (When the geometric inductance $L$ of the
SQUID is negligible, $\varphi_{+}=0$ and $\varphi_{-}=\pi$. As discussed in
Methods, a finite $L$ shifts $\varphi_{\pm}$ from these values Barone .)

We then have for the DC voltage

$V(\varphi_{\pm},I)=R_{\parallel}\sqrt{I^{2}-I_{\rm
cS}^{2}(\varphi_{\pm})},\quad(I>I_{\rm cS}(\varphi_{\pm}))$ (3)

with $I_{\rm cS}(\varphi_{\pm})=I_{0}\pm I_{\rm A}(V)$. As shown below,
allowing $I_{\rm A}$ and $R_{\parallel}$ to acquire a $V$ dependence yields a
close description of the measured curves $V(\varphi_{\pm},I)$.

Figure 1c shows $V(B,I)$ (measured in S1 at 135 mK) plotted in the $B$-$I$
plane. The black region ($V=0$) is bounded by the CPR curve $I_{\rm cS}(B)$.
At fixed $B$, $V$ is observed to increase steeply once $I$ excceeds $I_{\rm
cS}(B)$, approaching a linear increase at large $V$.

In Fig. 1d (main panel) we display a series of curves of $V$ vs. $I$ with $B$
as a parameter ($B\in[-24,\;-2]\;\mu$T). See Fig. S2a for curves with
$B\in[0,\;42]\;\mu$T. A series of steps are clearly seen in $V(B,I)$. As shown
in the inset, they persist to $V=100\,\mu$V. The steps lead to narrow peaks in
the derivative $dV/dI$ which follow a subharmonic sequence $V_{n}=V_{\rm
max}/n$. In Fig. 1e, we show that as $T$ is raised to 1.1 K, the peaks
converge to zero.

A key feature emerges when we investigate the effect of phase tuning on the
sequence of peaks. To see this, we plot the color map of $dV/dI(B,I)$ in the
$B$-$I$ plane (Fig. 2a). Above the CPR curve $I_{\rm cS}(B)$, we see a sheaf
of sharply defined sinusoidal curves that appear to peel off from the CPR
curve. Since $B$ is linear in the phase $\varphi$, we infer that each peak is
tracking a component of $I$ that varies sinusoidally with $\varphi$ – a
Josephson supercurrent.

We next show that each sinusoidal curve corresponds to an abrupt change in the
MAR order ($n\to n-1$) at a fixed $V_{n}$. Using the $I$-$V$ curves in Fig.
1d, we transform the vertical axis from $I$ to $V$. Under this transformation
(Fig. 2a$\to$2b), each sinusoidal curve in S1 collapses to a flat line. The
voltage $V_{n}$ on each line fits the subharmonic sequence $V_{n}=V_{\rm
max}/n$ for $n=2,3,4\cdots$ (Fig. S2c). The parameter $V_{\rm max}$ has the
$T$-independent value 200 $\mu$V below 300 mK but decreases as $T\to
T_{c}^{\rm Al}$ (1.2 K), consistent with $V_{\rm max}=2\Delta/e$ where
$\Delta$ is the energy gap of Al (Fig. S2b). These are key signatures of MAR
Klapwijk ; Averin ; Bratus ; Cuevas .

Increasing the junction spacing $d$ strongly damps both the critical current
$I_{\rm c}$ and the amplitude of the sinusoidal curves. In Fig. 3 we show
color maps of $dV/dI(B,I)$ for the four devices S1$\cdots$S4 measured at 340
mK. Although, at 340 mK, $I_{\rm c}$ (defined as half the peak-to-trough
excursion) decreases steeply with $d$ (see Fig. S2d in Methods), the
supercurrent is observable well beyond $d$ = 500 nm at 135 mK. Hereafter, we
focus on results from S1.

The striking staircase profile of $V(\varphi,I)$ (Fig. 1d) provides a vital
clue to the MAR charge transfer. When $V$ satisfies $V_{n}\leq V<V_{n-1}$,
where $V_{n}=2\Delta/ne$, the $n^{th}$-order process is dominant because of
the divergent DOS at the gap edge in Al. In an ideal junction ($\tau=1$), the
number of pairs shuttled is $n/2$ or $(n-1)/2$ for $n$ even or odd,
respectively (a quasiparticle is also transferred for odd $n$). In both cases,
the total charge shuttled is $ne$. Identifying each pair transfer as a
conductance channel, we have $n$ channels when $V_{n}\leq V<V_{n-1}$. At
either end of the interval, the channel number abruptly changes by 1.

Since $k_{B}T\ll\Delta/e$ at 135 mK, the normal current $I_{\rm N1}^{(n)}$ in
the shunt conductance, $G_{1}^{(n)}$, derives overwhelmingly from pairs that
are shuttled incoherently by MAR. As we show below (Fig. 4a), the shuttled
pairs give rise to both a normal current $I_{\rm N1}^{(n)}=VG_{1}^{(n)}$ and a
supercurrent; the former responds only to $V$ whereas the latter is sensitive
to $\varphi$.

From the staircase profile we infer a simple expression for $G_{1}^{(n)}$. For
$\tau=1$, the abrupt change in conductance at $V=V_{n}$, $\Delta
G_{1}=G_{1}^{(n+1)}-G_{1}^{(n)}$, equals a constant $g$ (the conductance for
one traversal). This implies $G_{1}^{(n)}=ng$.

In a non-ideal junction with $\tau<1$, classical scattering reduces the
transmission probability by the amount ($1-\tau$) at each reflection. Instead
of $G_{1}^{(n)}=ng$, we have

$G_{1}^{(n)}=g(1+\tau+\cdots+\tau^{n-1})=g(1-\tau^{n})/(1-\tau).$ (4)

As shown below, Eq. 4 leads to a quantitative description of
$V(\varphi_{\pm},I)$.

We turn next to the Josephson supercurrents $I_{s0}$ and $I_{s1}$. The
$n^{th}$-order sinusoidal curve in the $dV/dI$ map traces the variation of $I$
vs. $\varphi$ with $V=V_{n}$ (hereafter, we use $\varphi$ in place of $B$).
The maximum and minimum values attained by the current are called $I_{n+}$ and
$I_{n-}$, respectively. They occur at the extremal phases $\varphi_{+}$ and
$\varphi_{-}$, respectively.

By Eqs. 2 and 3, we have $I_{n\pm}^{2}=(V_{n}/R_{\parallel})^{2}+[I_{0}\pm
I_{\rm A}]^{2}$. Thus, the observed values $I_{n\pm}$ yield the observed
amplitudes $I_{\rm A}^{(n)}$, viz.

$I_{\rm A,obs}^{(n)}=[I_{n+}^{2}-I_{n-}^{2}]/(4I_{0}).$ (5)

These are plotted as black circles in Fig. 4a. Initially, $I_{\rm
A,obs}^{(n)}$ increases linearly with $n$ but curves downwards when $n>6$
(inset in Fig. 4a). To us, the initial $n$-linear growth of $I_{\rm
A,obs}^{(n)}$ is persuasive evidence that the supercurrent derives from the
population of Cooper pairs that are shuttled coherently.

At large $n$ ($V\to 0$), we expect $I_{\rm A,obs}^{(n)}$ to saturate since it
cannot exceed the critical current $I_{\rm c}$ at $V=0$. We propose that
saturation occurs because of the slight attenuation of the probability current
at each reflection, analogous to Eq. 4. Thus, $I_{\rm A}^{(n)}\sim
1+\tau+\cdots+\tau^{n-1}$ which we write as

$I_{\rm A}^{(n)}=I_{\rm c}[1-\tau^{n-1}].$ (6)

Using Eq. 6 to fit the data in Fig. 4a, we find that $I_{\rm c}$ = 27 $\mu$A,
close to $I_{\rm c}$ = 29.6 $\mu$A obtained in the CPR. The agreement supports
the reasoning behind Eq. 6. The value of $\tau$ is found to be 0.924, slightly
smaller than $\tau=$ 0.95 from the CPR fit (Eq. 11).

In Fig. 4a we also plot (as blue squares) the normal current in the sample
junction $I_{\rm N1}^{(n)}$, the product of $V_{n}$ and $G_{1}^{(n)}$ (for $g$
we use the value from the fit in Eq. 7). Despite the step-wise increase in
$G_{1}$ with $n$, the decrease of $V\sim 1/n$ forces $I_{\rm N1}^{(n)}$ to
decrease monotonically. By contrast, $I_{\rm A}^{(n)}$ (red and black symbols)
increases monotonically before saturating at $I_{\rm c}$. Both currents
reflect the series $1+\tau+\cdots+\tau^{n}$.

Heuristically, using Eqs. 4 and 6, we can adopt the generalized RSJ model to
describe the measured curves $V(\varphi_{\pm},I)$ by the expressions

$V(\varphi_{\pm},I)=\frac{\sqrt{I^{2}-[I_{0}\pm I_{\rm
A}^{(n)}]^{2}}}{G_{0}+G_{1}^{(n)}},$ (7)

where $n={\rm Int}[V_{\rm max}/V]$. The fits are shown as red curves in Fig.
4b. A more sensitive test is to compare the fits to the total observed
conductance $G_{\rm obs}(\varphi_{\pm},I)\equiv[I^{2}-(I_{0}\pm I_{\rm
A}^{(n)})^{2}]^{\frac{1}{2}}V(\varphi_{\pm},I)^{-1}$ (inset in Fig. 4b). As
seen, both curves of $G_{\rm obs}(\varphi_{\pm},I)$ fit well to Eq. 7. From
the fits, we find for $\varphi_{-}$ ($\varphi_{+}$)

$\tau=0.93\;(0.92),\quad G_{0}=3.00\;(2.94)\;{\rm S},\quad
g=0.32\;(0.20)\;{\rm S}.$ (8)

The value of $\tau$ is in agreement with the fit to Eq. 6. (For $\varphi_{+}$,
a discontinuous jump of $V$ occurs at the threshhold current $I_{\rm cS}$.
This implies a finite inertial term in the phase dynamics represented by a
shunt capacitor $C$, which is neglected for simplicity.)

Using the phase-tuning technique, we have uncovered a direct relation between
the voltage steps induced by MAR and the co-existing supercurrent. At large
$V$ (small $n$), the well-resolved steps in $V$ correspond to a step decrease
in the number of conductance channels for the pairs that are shuttled
incoherently. Even at $n=2$, the coherently shuttled pairs in Device S1
produce a supercurrent that is detectable (Fig. 2a). In the opposite limit
$V\to 0$, nearly all pairs are shuttled coherently. As a result, the amplitude
$I_{\rm A}$ saturates to $I_{\rm c}$ obtained in the CPR (Fig. 4a). In
between, the repeated reflections of the initial injected electron lays down
the MAR tracks. However, both the normal current observed as steps in $G_{1}$
and the co-existing supercurrent arise from the Cooper pairs that are shuttled
in its wake.

In high-transparency junctions that display MAR processes up to high order, we
expect a large fraction (possibly all) of the Josephson supercurrent to be
comprised of pairs coherently shuttled by the MAR. The ability to measure how
the amplitude $I_{\rm A}$ and other parameters ($n$, $\tau$, $g$, $\cdots$)
vary with $V$ in phase-tuned junctions may lead to a more quantitative
treatments of pairing correlations in the finite-$V$ regime, especially in
unconventional platforms SanJose ; vonOppen ; Halperin .

##  I Methods

###  I.1 Current definitions

We provide a glossary of the current parameters and definitions.
$I$ is the total current applied to the SQUID.
$I_{\rm cS}$ is the critical current of the SQUID.
$I_{0}$ is the critical current of the auxiliary junction.
$I_{\rm c}$ is the critical current of the sample junction when $V=0$.
$I_{\rm A}(V)$ is the $V$-dependent amplitude of the Josephson supercurrent at
finite $V$.
$I_{1}$ is the prefactor of the sample junction in Eq. 9.
$I_{n\pm}$ are the maximum and minimum values of $I$ in the sinusoidal curve
with $V=V_{n}$.
$I_{\rm N1}$ and $I_{\rm N0}$ are the normal currents in the sample and
auxiliary junctions, respectively, at finite $V$.
$I_{\rm N}$ is the total normal current in the SQUID.
$I_{s0}$ is the supercurrent in the auxiliary junction.
$I_{s1}$ is the supercurrent in the sample junction.

###  I.2 Device fabrication and measurement

We used double-layer e-beam lithography, in combination with tilted-substrate
thermal evaporation, to fabricate the $S$-${\cal I}$-$S$ junction. Initially,
the substrate is spin-coated with MMA EL11 at 3 000 rpm for 30 s twice and
baked at 175∘ C for 5 min, followed by spin-coating with PMMA 950 A07 at 4 000
rpm for 60 s and then baked at 175∘ C for 5 min. Next, the SQUID pattern was
e-beam written using the Raith eLiNE writer with beam energy set at 30 kV,
aperture at 10 $\mu$m and the dose level at $\sim$300 $\mu$C/cm2. After
developing in MIBK solution (MIBK: IPA=1:3) for 3 min and rinsing in IPA
solution for 1 min, we fabricated a suspended bridge, using a half-dose beam
to remove the underlying MMA layer, while keeping the suspended upper PMMA
layer intact. The chip with the pattern defined is placed inside a thermal
evaporator equipped with a tiltable stage, with vacuum at $\sim 5\times
10^{-7}$ mbar. To remove residual resist and several (oxidized) monolayers of
MoTe2, we exposed the chip to an RF Argon plasma _in situ_. After cleansing,
the first layer of Al (60 nm) is deposited at a rate of 10 $\AA$/s. Then a
mixture of Ar/O2 (10$\%$ O2) is injected into the chamber for 30 min at
$5\times 10^{-3}$ mbar to oxidize the Al. A second layer of Al (120 nm) is
next deposited with the angle set at a new value to define an overlapping Al-
AlOx-Al junction under the suspended PMMA resist bridge. Finally, the device
is immersed in acetone to wash off the extra Al layer.

The $I$-$V$ measurements were performed in a top-loading wet dilution
refrigerator (Oxford Instruments Kelvinox TLM400) with base temperature of
15–20 mK. The dc bias and ac excitation current is provided by Agilent 33220A
function generator with a bias resistance of 100 k$\Omega$. After
preamplification (NF LI-75A preamp.), the SQUID voltage $V$ was fed to a lock-
in amplifier (Stanford Research SR830) for $dV/dI$ measurement, as well as a
nanovoltmeter (Keithley 2182A) for $I-V$ measurement. Data above 300 mK were
acquired in a Heliox Helium-3 cryostat with base temperature of 340 mK.

###  I.3 The CPR curve

We describe the fitting procedure to obtain the current-phase relation (CPR)
curve from measurements on Device S1 at 135 mK. The CPR is the curve of the
critical current $I_{\rm cS}(\varphi)$ bounding the dissipationless region
where $V=0$. The total supercurrent is written as

$I(B,\delta)=I_{0}\sin(\delta+\varphi)+I_{1}\frac{\sin\delta}{\sqrt{1-\tau\sin^{2}(\delta/2)}},$
(9)

with $\varphi=2\pi\phi/\phi_{0}$ where $\phi$ is the total flux piercing the
SQUID and $\phi_{0}=h/2e$ is the superconducting flux quantum. For the $SNS$
junction, we adopted the CPR expression ${\cal
F}(\delta)\equiv\sin\delta/\sqrt{1-\tau\sin^{2}(\delta/2)}$ appropriate for
subgap states in junctions with high transparency ($\tau\to 1$). Its prefactor
$I_{1}$ is proportional to $I_{\rm c}$. The total flux $\phi$ is the sum of
the applied flux $\phi_{\rm a}$ ($\phi_{\rm a}=BA$) and the Amperean flux
“$LI$” produced by supercurrents flowing around the SQUID loop. We have, when
$V=0$,

$\phi=\phi_{a}-L_{0}I_{0}\sin(\delta+\varphi)+L_{1}I_{1}{\cal F(\delta)}.$
(10)

The partial equivalent inductances $L_{0}$ and $L_{1}$ represent the
partitioning of $L$ into contributions from the two branches ($L=L_{1}+L_{2}$)
Fulton .

In our device layout, we have $L_{0}\simeq L_{1}\simeq L/2$. We estimate the
total inductance $L$ as that of a square loop (of sides 6 $\mu$m and linewidth
1 $\mu$m) to be 8.2 pH. (This is close to the value estimated from the shift
of the mimima of each sinusoidal curve in the color map of $dV/dI$ (Fig. 2a)
$dI/dB$ which yields $L\sim$ 9 pH as shown below.)

With the total flux given by Eq. 10, the optimization of the fit by analytical
means is highly unstable. Instead, we used the following numerical procedure.
Starting with seed values for the 3 unknowns, $\tau$, $I_{0}$ and $I_{1}$, we
allow $\varphi$ and $\delta$ to range over the full parameter space
$\varphi\in[0,2\pi]$, $\delta\in[0,2\pi]$. The magnitude of
$I(\varphi,\delta)$ at each point $(\varphi,\delta)$ (computed using Eqs. 9
and 10 with $L$ = 8.2 pH) defines a surface $I(\varphi,\delta)$. When we
project the surface onto the $\varphi$-$I$ plane, its upper boundary (maximum
$I$) gives the CPR curve (using Eq. 10 to convert $\varphi$ to $B$). The
deviation of this calculated CPR from the observed CPR yields an error
function (a scalar map in the space of $\tau,I_{0},I_{1}$). By iteration, we
converge rapidly to the optimal fit values.

$\tau=0.95,\quad I_{0}=93\,\mu{\rm A},\quad I_{1}=24\,\mu{\rm A}.$ (11)

In the $V=0$ regime, we define the critical current $I_{\rm c}$ of the
$S$-$N$-$S$ junction as one-half the trough-to-peak excursion of the CPR
curve. From the red circles in Fig. S1, we get $I_{\rm c}\simeq$ 29.6 $\mu$A.
Hence the prefactor $I_{1}\simeq 0.81\;I_{\rm c}$.

###  I.4 Data in Device S1

Figure S2a shows curves of $V(B,I)$ in Device S1 for $0<B<42\;\mu$T
(complements the curves in Fig. 1a). Figure S2b shows that $V_{\rm gap}$ in S1
has a $T$ dependence consistent with the gap order parameter in Al. In S2c,
the plot verifies that $V_{n}$ measured at the flat lines in Fig. 2b satisfies
$V_{n}=V_{\rm max}/n$. Both Panels (b) and (c) support the identification of
the steps in $V(B,I)$ with MAR. In Fig. S2d, we show the decrease in $I_{\rm
c}$ with junction spacing $d$ ($I_{\rm c}$ is inferred from the CPR curves).
At $T$= 135 mK, the decrease is gradual, but at 340 mK, $I_{\rm c}$ decays
steeply for $d>$200 nm.

###  I.5 RSJ model of single junction

First we apply the RSJ model to the $SNS$ junction 1 in isolation. Treating it
as an overdamped Josephson junction (JJ) (setting the inertial term $C=0$), we
have a resistance $R$ in parallel with the JJ. When the applied current
$I>I_{\rm c}$, winding of the phase $\theta$ leads to a time-dependent voltage
$V(t)=\hbar\dot{\theta}(t)/2e$. Adding the Josephson supercurrent to
$V(t)/R_{1}$, the total current $I$ is

$I=\hbar\dot{\theta}/2eR+I_{c}\sin\theta.$ (12)

The differential equation may be integrated to solve for the winding rate

$\dot{\theta}(t)=\frac{p^{2}-1}{p}\left(\frac{\cos\theta}{\cos\frac{\tau}{2}\sqrt{p^{2}-1}}\right)^{2}$
(13)

where $\tau=(2I_{c}eR/\hbar)t$, and $p=I/I_{c}$. The winding rate
$\dot{\theta}(t)$ is a periodic function of $t$ with narrow peaks separated by
the period $T$ given by

$T=\frac{2\pi}{\sqrt{I^{2}-I_{\rm c}^{2}}}\frac{\hbar}{2eR}.$ (14)

Time-averaging the winding rate $\langle\dot{\theta}\rangle$, we obtain for
the DC voltage

$\langle
V\rangle=\frac{\hbar}{2e}\oint\dot{\theta}\frac{dt}{T}=R\sqrt{I^{2}-I_{\rm
c}^{2}}$ (15)

###  I.6 Asymmetric SQUID with finite inductance $L$

In the finite-$V$ regime, we have for the total current

$I=\frac{\hbar\dot{\gamma}}{2eR_{\parallel}}+I_{0}\sin\gamma+I_{\rm
A}(V)\sin\delta,$ (16)

where we assume that ${\cal F}(\delta)\simeq\sin\delta$ and replace $I_{\rm
c}$ with the amplitude $I_{\rm A}(V)$. The color plot in Fig. 2a shows that
the sinusoidal assumption is reasonable except very close to the CPR curve
($V\to 0$). In applied $B$, the phases $\gamma$ and $\delta$ remain related by
the total flux $\phi$:

$\gamma-\delta=\varphi=\frac{2\pi}{\phi_{0}}\phi.$ (17)

In the presence of the normal currents $I_{\rm N0}=V/R_{0}$ and $I_{\rm
N1}=V/R_{1}$ driven by $V$, Eq. 10 changes to

$\phi=\phi_{a}-L_{0}I_{0}(\sin\gamma+\frac{I_{\rm N0}}{I_{0}})+L_{1}I_{\rm
A}(\sin\delta+\frac{I_{\rm N1}}{I_{\rm A}}).$ (18)

Defining $\varphi_{a}=(2\pi/\phi_{0})\phi_{a}$ and eliminating $\phi$ using
Eq. 17, Eq. 10 becomes

$\varphi_{a}=\gamma-\delta+\beta_{0}(\sin\gamma+\frac{I_{\rm
N0}}{I_{0}})-\beta_{1}(\sin\delta+\frac{I_{\rm N1}}{I_{\rm A}}).$ (19)

with the beta parameters $\beta_{0}=2\pi I_{0}L_{0}/\phi_{0}$ and
$\beta_{1}=2\pi I_{\rm A}L_{1}/\phi_{0}$. In Device 1, $\beta_{0}\simeq 1.17$
and $\beta_{1}\simeq 0.38$.

Close to the threshold ($V\to 0$), the nonlinear equations Eqs. 16 and 19
require a numerical solution, but the solutions are difficult to relate to
measurable quantities. We restrict attention to the voltage regime $V>V_{\rm
max}/n_{\rm max}$ where $n_{\rm max}=11$ is the highest MAR order resolved. In
this limit, $I_{\rm N1}\simeq 0.9I$ and $I_{\rm N0}\simeq 0.1I$, i.e. $90\%$
of the applied current flows through the SNS junction as a normal current.
Then Eq. 19 simplifies to

$\varphi_{a}\to\gamma-\delta+\beta_{0}\sin\gamma-\beta_{1}\frac{I_{\rm
N1}}{I_{\rm A}}.$ (20)

Substituting this into Eq. 16, we see that at $\varphi_{+}$, when the
supercurrents are parallel ($I_{s0}\parallel I_{s1}$ ) the corresponding
$\varphi_{a}$ is

$\varphi_{a+}=\beta_{0}\sin\gamma-\beta_{1}\frac{I_{\rm N1}}{I_{\rm A}}.$ (21)

This shows that the maximum $I_{+}^{(n)}$ in a sinusoidal curve in Fig. 2a is
shifted from 0 by an amount linear in $I_{\rm N1}=0.9I$. Similarly,
$\varphi_{-}$ is shifted from $\pi$ by the same amount as in Eq. 21.

In the large $I$ limit, Eq. 21 further simplifies to

$\varphi_{a+}\simeq\frac{2\pi}{\phi_{0}}\frac{LV}{2R_{1}}.$ (22)

The $V$-linear shifts of $\varphi_{\pm}$ are clearly observed in the color map
of $dV/dI$ (Fig. 2a).

The primary effect of $L$ is to induce a phase shift in the sinusoidal curves
without affecting the amplitude. In our fits to Eq. 7, $L$ does not enter
explicitly because we use the experimentally observed shifts to locate
$\varphi_{a+}$ and $\varphi_{a-}$, respectively.

## References

  * (1) A. F. Andreev, The Thermal Conductivity of the Intermediate State in Superconductors, _Zh. Eksp. Teor. Fiz._ 46, 1823 (1964) [Sov. Phys. JETP 19, 1228 (1964)].
  * (2) Hubert B. Heersche, Pablo Jarillo-Herrero, Jeroen B. Oostinga, Lieven M. K. Vandersypen and Alberto F. Morpurgo, Bipolar supercurrent in graphene, _Nature_ 446, 56-59 (2007). doi:10.1038/nature0555
  * (3) Xu Du, Ivan Skachko, and Eva Y. Andrei, Josephson current and multiple Andreev reflections in graphene SNS junctions, Phys. Rev. B77, 184507 (2008). DOI: 10.1103/PhysRevB.77.184507
  * (4) J-D. Pillet, C. H. L. Quay, P. Morfin, C. Bena, A. Levy Yeyati and P. Joyez, Andreev bound states in supercurrent-carrying carbon nanotubes revealed, _Nature Physics_ 6, 965-969 (2010). DOI: 10.1038/NPHYS1811
  * (5) M. Zgirski, L. Bretheau, Q. Le Masne, H. Pothier, D. Esteve, and C. Urbina, Evidence for Long-Lived Quasiparticles Trapped in Superconducting Point Contacts, Phys. Rev. Lett. 106, 257003 (2011). DOI: 10.1103/PhysRevLett.106.257003
  * (6) L. Bretheau, C. O. Girit, C. Urbina, D. Esteve, and H. Pothier Supercurrent Spectroscopy of Andreev States, _Phys. Rev. X_ 3, 041034 (2013) DOI: 10.1103/PhysRevX.3.041034
  * (7) E. Strambini, S. D’Ambrosio, F. Vischi, F. S. Bergeret, Yu. V. Nazarov and F. Giazotto, The $\omega$-SQUIPT as a tool to phase-engineer Josephson topological materials, _Nature Nanotechnology_ 11, 1055–1059 (2016). DOI: 10.1038/NNANO.2016.157
  * (8) Anil Murani _et al._ , Ballistic edge states in Bismuth nanowires revealed by SQUID interferometry, _Nat. Commun._ 8, 15941 (2017). doi: 10.1038/ncomms15941.
  * (9) Sean Hart, Hechen Ren, TimoWagner, Philipp Leubner, Mathias Mühlbauer, Christoph Brüne, Hartmut Buhmann, Laurens W. Molenkamp and Amir Yacoby, Induced superconductivity in the quantum spin Hall edge, _Nat. Phys._ 10 638-643 (2014). DOI: 10.1038/NPHYS3036
  * (10) Wudi Wang, Stephan Kim, Minhao Liu, F. A. Cevallos, R. J. Cava, N. P. Ong, Evidence for an edge supercurrent in the Weyl superconductor MoTe2, _Science_ 368, 534–537 (2020). 10.1126/science.aaw9270.
  * (11) A. Furusaki, H. Takayanagi and M. Tsukada, Theory of quantum conduction of supercurrent through a constriction, Phys. Rev. Lett. 67, 132-135 (1991).
  * (12) C. W. J. Beenakker and H. van Houten, Josephson current through a superconducting quantum point contact shorter than the coherence length, Phys. Rev. Lett. 66, 3056-3059 (1991).
  * (13) Akira Furusaki, Josephson current carried by Andreev levels in superconducting quantum point contacts, _Superlattice and Microstructures_ 25, 809-818 (1999).
  * (14) Pablo San-Jose, Jorge Cayao, Elsa Prada and Ramon Aguado, Multiple Andreev reflection and critical current in topological superconducting nanowire junctions, _New Jnl. Phys._ 15, 075019 (2013). doi:10.1088/1367-2630/15/7/075019
  * (15) Yang Peng, Falko Pientka, Erez Berg, Yuval Oreg, and Felix von Oppen, Signatures of topological Josephson junctions, _Phys. Rev._ B 94, 085409 (2016). DOI: 10.1103/PhysRevB.94.085409
  * (16) Falko Pientka, Anna Keselman, Erez Berg, Amir Yacoby, Ady Stern, and Bertrand I. Halperin, Topological Superconductivity in a Planar Josephson Junction, _Phys. Rev. X_ 7, 021032 (2017). DOI: 10.1103/PhysRevX.7.021032
  * (17) B. N. Taylor and E. Burnstein, Excess currents in electron tunneling between superconductors, Phys. Rev. Lett. 10, 14 (1963).
  * (18) S. M. Marcus, The magnetic field dependence of the 2$\Delta/n$ structure observed in Pb-PbO-Pb superconducting tunneling junctions, _Phys. Lett._ 19, 623 (1966).
  * (19) T. M. Klapwijk, G. E. Blonder and M. Tinkham, Explanation of subharmonic energy gap structure in superconducting contacts, _Physica_ 109-110B, 1657-1664 (1982).
  * (20) D. Averin and A. Bardas, ac Josephson Effect in a single quantum channel, Phys. Rev. Lett. 75, 1831-1834 (1995).
  * (21) E. N. Bratus, V. S. Shumeiko, and G. Wendin, Theory of subharmonic gap structure in superconducting mesoscopic tunnel contacts, Phys. Rev. Lett. 74, 2111-2113 (1995).
  * (22) J. C. Cuevas, A. Martin-Rodero and A. Levy Yeyati, Hamiltonian approach to the transport properties of superconducting quantum point contacts, Phys. Rev. B54, 7366-7379 (1996).
  * (23) A. Th. A. M. De Waele and R. De Bruyn Ouboter, Quantum intereference phenomena in point contacts between two superconductors, _Physica_ 41, 225 (1969).
  * (24) T. A. Fulton, L. N. Dunkleberger, and R. C. Dynes Quantum interference properties of double Josephson Junctions, Phys. Rev. B6, 855 (1972).
  * (25) A. Barone and G. Paterno, _Applications and Physics of the Josephson Effect_ (John Wiley and Sons, 1982), Ch. 12, p. 375.

* * *

§Corresponding author email<EMAIL_ADDRESS>

Acknowledgements
S. K. and the experiments performed at milliKelvin temperatures were supported
by an award from the U.S. Department of Energy (DE-SC0017863). The crystal
growth, led by R. J. C., L. M. S. and S. L., was supported by MRSEC award from
the U.S. National Science Foundation (NSF DMR-2011750). Z. Y. Z. and N. P. O.
were supported by the Gordon and Betty Moore Foundation’s EPiQS initiative
through grant GBMF4539.

Author contributions
The experiment was designed by Z.Y.Z. and N.P.O. and carried out by Z.Y.Z. and
S.K. on successive generations of MoTe2 crystals grown by S.L., L.M.S. and
R.J.C. Analysis of the data was done by Z.Y.Z. and N.P.O. The manuscript was
written by N.P.O. and Z.Y.Z.

Additional Information
Supplementary information is available in the online version of the paper.
Correspondence and requests for materials should be addressed to N.P.O.

Competing financial interests
The authors declare no competing financial interests.

Figure Captions
Figure 1
Multiple Andreev reflections in the asymmetric SQUID layout with the sample
$S$-$N$-$S$ junction based on MoTe3. Panel (a) shows a sketch of the
$5^{th}$-order MAR process. Right-moving electrons (red circles) are Andreev
reflected as left-moving holes (white circles). After 5 traversals, the
excitation gains sufficient energy to scale the gap barrier. In the process 2
pairs (green arrows) are shuttled across. The density of states (DOS) in Al
are shaded grey. Panel (b) shows four devices fabricated on an exfoliated
crystal of MoTe2 (gray area). The sketch in the upper inset shows the
$S$-$N$-$S$ ($S$-${\cal I}$-$S$) junction on the right (top) branch. The
applied flux $\phi_{\rm a}=BA$ pierces the enclosed rectangle (see scale bar).
In Panel (c) the color map is comprised of 150 $I$-$V$ curves with spacing
$\Delta B=2\,\mu$B (color scale in vertical bar) measured in S1 at $T$ = 135
mK. At fixed $B$, $V$ increases steeply once $I$ exceeds the SQUID critical
current $I_{\rm cS}(B)$ (the CPR curve). If $V$ is held fixed, the current
varies periodically with vs. $B$ up to $V\sim$ 80 $\mu$V (thin black curves
are contours of $V(B,I)$). Panel (d) (main panel) shows the $I$-$V$ curves
with $B$ varying from -24 $\mu$T to -2 $\mu$T in steps of 2 $\mu$T.
Subharmonic steps ($n$ as indicated) are visible up to $V$ = 105 $\mu$V (see
inset). Peaks in $dV/dI$ at the steps occur at $V_{n}=V_{\rm max}/n$
($n=2,3,\cdots$). Panel (d) shows traces of $dV/dI$ vs. $V$ for $0.34<T<1.10$
K. As $T\to T_{\rm c}^{\rm Al}$, $V_{\rm max}\to 0$ in accord with MAR.

Figure 2
Color maps of the differential resistance $dV/dI(B,I)$ measured in S1 at 135
mK, displayed in the $(B,I)$ plane (Panel a) and the $(B,V)$ plane (b) with
scale bar on right. In Panel a, the black area ($dV/dI=0$) is bounded by the
CPR curve $I_{\rm cS}(B)$. Above the CPR, the series of sinusoidal curves
trace the variation of narrow peaks in $dV/dI$ vs. $B$. In Panel b, the
vertical axis is transformed to $V$ using the $I$-$V$ curve at each $B$. Each
sinusoidal curve collapses to a flat line with $V$ fixed at the MAR
$V_{n}=V_{\rm max}/n$.

Figure 3
Comparison of color maps of $dV/dI(I,B)$ measured in S1–S4 at $T$ = 340 mK.
The junction spacing $d$ is 200, 300, 400, 500 nm in S1 to S4, respectively.
As $d$ increases, the supercurrent amplitude $I_{\rm c}$ inferred from the CPR
decreases rapidly at 340 mK (see Fig. S2d in Methods). Unlike at 135 mK,
$dV/dI$ is not strictly zero below the CPR curve because of thermally excited
quasiparticles.

Figure 4
Fits to the supercurrent amplitude $I_{\rm A,obs}^{(n)}$, the $I$-$V$ curves
and total conductance $G_{\rm obs}$ in Device S1. Panel (a): The observed
supercurrent amplitude $I_{\rm A,obs}^{(n)}$ vs. $n$ (Eq. 5, black circles)
compared with $I_{\rm A}^{(n)}$ calculated using Eq. 6 (red circles). The best
fit is obtained with $I_{\rm c}$ = 27 $\mu$A and $\tau=0.924$. The normal
current $I_{\rm N1}^{(n)}=V_{n}G_{1}^{(n)}$ is plotted as blue squares. The
data for $I_{\rm A,obs}^{(n)}$ and fit are shown in expanded view in the
inset. In Panel (b) (main panel), we compare the measured $V(\varphi_{\pm},I)$
(black curves) with the best fit to Eq. 7 (red curves). The integer $n$ at
each step is indicated. For clarity, the curve $V(\varphi_{+},I)$ is shifted
upwards by 12.5 $\mu$A. The inset shows the same fits compared against the
observed conductances $G_{\rm obs}(\varphi_{\pm},I)$ defined below Eq. 7,
which provide a more sensitive test for the fits.

Figure S1
Numerical procedure for fitting the CPR equation Eq. 9 to the curve of $I_{\rm
cS}(B)$ measured in S1 at 135 mK. The $x$ axis is first converted to the
corrected field $B_{\rm cor}$ using Eq. 10. The set of values for $I(B_{\rm
cor},\delta)$ is calculated from Eq. 9 allowing $\delta$ and $\varphi$ to vary
over $[0,2\pi]\times[0,2\pi]$. These values are projected onto the $I$-$B_{\rm
cor}$ plane. The difference between the maximum value taken by $I$ and the
measured CPR (red circles) defines the error function. We iteratively reduce
the error function to force convergence to the optimal values Eq. 11. The wavy
motif in gray is an artefact of the procedure with no physical significance.

Figure S2
Panel (a) displays $I$-$V$ curves with $B\in[0,42]\;\mu$T in steps of $\Delta
B=2\;\mu$T (complement of Fig. 1d). The integers $n$ at each step is
indicated. Panel (b) plots the $T$ dependence of $V_{\rm gap}$ inferred from
the MAR sequence at each $T$. The $T$ dependence is consistent with the gap
profile in Al. Panel (c) shows the plot of $V_{n}$ vs. $1/n$. The linear
relationship $V_{n}=V_{\rm max}/n$ confirms the identification of the steps
with MAR. Panel d: The variation of the $S$-$N$-$S$ critical current $I_{\rm
c}$ vs. junction spacing $d$ at 135 mK (red symbols) and 340 mK (black). At
135 mK, $I_{\rm c}$ asymptotes to $\sim 30\;\mu$A as $d$ decreases below 200
nm.

Figure 1:  Figure 2:  Figure 3:  Figure 4:  Figure S1:  Figure S2:


Confining type-II spherical core-shell quantum dot heterostructures with
narrow and wide band gaps

T<EMAIL_ADDRESS>M.S.
Nurisya†††risya@upm.edu.mya,b and A. Jellal‡‡‡a.jellal@ucd.ac.mac,d

aDepartment of Physics, Faculty of Science, Universiti Putra Malaysia,
43400 UPM Serdang, Selangor, Malaysia

bLaboratory of Computational Sciences & Mathematical Physics,
Institute for Mathematical Research, Universiti Putra Malaysia, 43400 UPM
Serdang, Selangor, Malaysian

cLaboratory of Theoretical Physics, Faculty of Sciences, Chouaïb Doukkali
University,

PO Box 20, 24000 El Jadida, Morocco

dCanadian Quantum Research Center, 204-3002 32 Ave Vernon,
BC V1T 2L7, Canada

Using a single-band model, the lowest transition energy was analysed between
the lowest unoccupied molecular orbital (LUMO) of the conduction band and the
highest occupied molecular orbital (HOMO) of the valence band. We focus on
categorising the confinement strength in type-II core-shell quantum dots
(CSQDs) based on the step-potential and show how it will affect their
transition energy. Our model is applied to CSQDs of the heterostructures
PbS/CdS and ZnTe/ZnSe through narrow and wide band gaps, respectively. It
found that PbS/CdS CSQDs demonstrates a strong confinement in which their
transition energy would increase more compared to its weak confinement case in
ZnTe/ZnSe CSQDs. The weak confinement case also demonstrated both blue-shift
and red-shift of photoluminescence emission compared to the bulk ZnTe and ZnSe
for which it can be inferred as pseudo type-II CSQDs. This would help
experimentalists to tune the transition energy of type-II model in order to
fabricate photons with longer carrier lifetime compared to the type-I model.

Keywords: Core-shell quantum dots, transition energy, step-potential, low
dimensional semiconductor

##  1 Introduction

Quantum dots realize the limiting case of size quantization in semiconductors,
which cause modifications of their electronic properties [1]. The colloidal
core-shell quantum dots (CSQDs) are heterostructures and can be obtained by
chemical synthesis [2, 3, 4] with reproducible and controllable size and shape
and low fabrication cost. The demand and applications of nanotechnology in
devices for biological tagging [5, 6, 7, 8], diodes [7], solar cells [5, 7,
9], and lasers [7, 8, 9] have increased over the years made quantum dots to be
studied extensively, and developed further into core-shell structures. Through
this chemical surface passivation, core-shell quantum dots (CSQDs) would have
better photo-stability [10] and a highly tunable optical spectra [10, 11].
Capping a quantum dots surface with another semiconductor material would also
give room for wavefunction engineering. Then it is possible to form CSQDs with
type-I band alignments or type-II band alignments for the same pair of
materials by manipulating the materials composition or by altering their
environment by introducing electric field [6, 12].

While a type-I CSQDs would have both of their charge carriers confined within
the same region, in type-II CSQDs, both electrons and holes would be confined
in different regions; one within the core region and the other within the
barrier interface region [13]. This novel ability to control the region of
carriers confinement leads to lower non-radiative recombination rate [14] in
type-II CSQDs. Type-II also offers longer carrier lifetimes while reducing the
radiative recombination probability [15, 16]. In addition, type-II CSQDs also
preferred in engineering a near-infrared or infrared photons that can be used
to assist surgery for direct visual guidance by injecting the photons into the
skin [17].

Due to these promising features of type-II CSQDs, it is best to acknowledge
the importance to understand the effects of step-potentials in type-II CSQDs
by using the single band model and see how it will affect their transition
energy. This is much due to the previous findings on how different confinement
strength in type-I CSQDs affects their transition energy [18]. Previous works
have tackled the transition energy problem in type-II CSQDs using continuum
mechanical model (CMM), valence force field model (VFF) and effective mass
approximation (EMA) [19, 20, 21]. However, CMM and VFF approaches involved
8-band k$\cdot$p method while EMA often involved with (2,8)-band Hamiltonian
approach which often hard to solve and require longer calculations compared to
the single-band model [19, 21]. As the name suggests, the single-band model
only focusing on the transition energy between the lowest unoccupied molecular
orbital (LUMO) of conduction band and the highest occupied molecular orbital
(HOMO) of valence band.

The present work reports the computation of the lowest transition energy
between the lowest unoccupied molecular orbital of the conduction band and the
highest occupied molecular orbital of the valence band. This will be done by
using a single band model involving the confining potentials for electrons and
holes. The solutions of the energy spectrum will be determined by requiring
separable wave functions in radial and angular parts. The radial part will
involve the physical parameters generated from the confinements.

The paper is organized as follows. Section 2 would be dedicated for discussion
on the theoretical single-band approach on energy transition observed in
CSQDs. In section 3, the mathematical formulation of the single-band theory is
applied on type-II CSQDs for strong and weak confinement cases for which
analogous to the CSQD heterostructure with narrow and wide band gap. The
analysis of their radial probability will be presented in section 4. We close
by giving conclusion remarks.

##  2 Theory

We consider type-II core-shell quantum dots (CSQDs) based on the step-
potential. Figure 1 shows a schematic diagram for type-II CSQDs, with
$E_{c(v)}$ and $V_{c(v)}$ are the energy and step-potential for
conduction(valence) band, while $r_{o}$ and $r_{s}$ are the radius of core
quantum dots and CSQDs respectively. For carriers outside the shell layer, the
corresponding $V(r)$ would be infinite. This assumption of having an infinite
barrier would lead to a complete spatial separation between the CSQDs and its
environment. Therefore the wavefunction overlap and the Coulomb potential
outside the “artificial atom” are negligible [13].

Figure 1: (color online) Schematic diagram of band line-up of Type-II CSQDs.

To achieve our goal, let us introduce the step-potential due to electrons
confinement [18]

$\displaystyle V_{e}(r)=\left\\{\begin{matrix}0&;&0\leq r\leq r_{o}\\\
V_{c}&;&r_{o}\leq r\leq r_{s}\\\ \infty&;&r\geq r_{s}.\end{matrix}\right.$ (1)

and focus only interesting solution comes from solving the radial part
$R_{nl}(r)$ of the following three-dimensional wavefunction

$\Psi_{nlm}(\vec{r})=R_{nl}(r)Y_{lm}(\theta,\phi),$ (2)

where $Y_{lm}$ are its spherical harmonic functions, with $n$, $l$ and $m$ are
the principal, orbital and magnetic quantum numbers respectively. Since this
work is focusing on analysing a single-band model, hence one only needs to
consider the transition energy between the lowest unoccupied molecular orbital
(LUMO) of the conduction band and the highest occupied molecular orbital
(HOMO) of valence band, that is taking $n=1$ and $l=m=0$. According the
potential profile Eq. (1), the radial wavefunction for ground state electron
then would be

$\displaystyle R_{i}(r)=\left\\{\begin{matrix}R_{1}(r)=AJ_{l}(k_{1}r)&;&0\leq
r\leq r_{o}\\\ R_{2}(r)=BH_{l}(k_{2}r)+CN_{l}(k_{2}r)&;&r_{o}\leq r\leq
r_{s}\\\ R_{3}(r)=0&;&r\geq r_{s}.\end{matrix}\right.$ (3)

where $j_{l}(k_{i})$, $h_{l}(k_{i})$ and $n_{l}(k_{i})$ are spherical Bessel
function of the first kind, modified spherical Bessel function of the first
kind and modified spherical Bessel function of the second kind respectively,
with $i=1$ for core, $i=2$ for shell and $i=3$ for outside shell layer. The
involved wave vector for both regions are

$k_{1}=\frac{\sqrt{2m^{e}_{1}E_{e}}}{\hbar},\qquad
k_{2}=\frac{\sqrt{2m^{e}_{2}[V_{c}(r)-E_{e}]}}{\hbar}$ (4)

and $m_{i}^{e}$ is the effective mass of electrons. Eq. (3) would form a
transcendental equation that obey det$(E_{e})=0$ and thus

$\dfrac{m_{2}^{*}}{m_{1}^{*}}\dfrac{J^{\prime}_{0}(k_{1}r_{o})}{J_{0}(k_{1}r_{o})}=\dfrac{H^{\prime}_{0}(k_{2}r_{o})n_{0}(k_{2}r_{s})-H_{0}(k_{2}r_{s})N^{\prime}_{0}(k_{2}r_{o})}{H_{0}(k_{2}r_{o})n_{0}(k_{2}r_{s})-H_{0}(k_{2}r_{s})N_{0}(k_{2}r_{o})}.$
(5)

One will obtain the energy of electrons, $E_{e}$ after solving Eq. (5). As for
these electron states, the numerical analysis is the same as solving for
electron states in Type-I CSQDs, as their electrons are confined within the
conduction band of the core region.

In the case of Type-II, the confinement for holes would be no longer within
the valence band of core region, but instead, it would be spatially confined
within the valence band of both core and shell regions. This would reflect
different step-potential distribution compared to their electron counterpart,
where they only confined within their core region as

$\displaystyle V_{h}(r)=\left\\{\begin{matrix}V_{v}&;&0\leq r\leq r_{o}\\\
0&;&r_{o}\leq r\leq r_{s}\\\ \infty&;&r\geq r_{s}\end{matrix}\right.$

and then the corresponding radial wavefunctions for ground state of holes are

$\displaystyle
R_{i}(r)=\left\\{\begin{matrix}R_{1}(r)=DT_{l}(k^{\prime}_{1}r)&;&0\leq r\leq
r_{o}\\\ R_{2}(r)=FP_{l}(k^{\prime}_{2}r)+GQ_{l}(k^{\prime}_{2}r))&;&r_{o}\leq
r\leq r_{s}\\\ R_{3}(r)=0&;&r\geq r_{s},\end{matrix}\right.$ (6)

where $t_{l}(k_{i})$, $p_{l}(k_{i})$ and $q_{l}(k_{i})$ are modified spherical
Bessel function of the first kind, spherical Bessel function of the first kind
and spherical Bessel function of the second kind, respectively. The wave
vectors are

$k^{\prime}_{1}=\frac{\sqrt{2m^{h}_{1}[V_{c}(r)-E_{h}]}}{\hbar},\qquad
k^{\prime}_{2}=\frac{\sqrt{2m^{h}_{2}E_{h}}}{\hbar}$ (7)

and $m_{i}^{h}$ is the effective mass of holes. Also here Eq. (6) will form
the transcendental equations that can be written as

$\dfrac{m_{2}^{h}}{m_{1}^{h}}\dfrac{T_{0}(k^{\prime}_{1}r_{o})}{T_{0}(k^{\prime}_{1}r_{o})}=\dfrac{P^{\prime}_{0}(k^{\prime}_{2}r_{o})Q_{0}(k^{\prime}_{2}r_{s})-P_{0}(k^{\prime}_{2}r_{s})Q^{\prime}_{0}(k^{\prime}_{2}r_{o})}{P_{0}(k^{\prime}_{2}r_{o})Q_{0}(k^{\prime}_{2}r_{s})-P_{0}(k^{\prime}_{2}r_{s})Q_{0}(k^{\prime}_{2}r_{o})},$
(8)

which would give the numerical value for energy of the holes, $E_{h}$.

Combining all solutions, the lowest transition energy, $E_{\text{qd}}$ for
CSQDs can be solved by using the Brus equation [22]

$E_{\text{qd}}=E_{g}^{bulk}+E_{e}+E_{h}-\frac{1.8e^{2}}{4\pi\varepsilon_{o}\varepsilon_{r}r_{o}},$
(9)

where $e$ is the elementary charge, $\epsilon_{o}$ is the permittivity of free
space, $\epsilon_{r}$ is the relative permittivity and $r_{o}$ is the core
radius.

##  3 Results and discussions

###  3.1 Type-II CSQDs with strong confinement

Now the obtained mathematical formulation will be applied to a type-II CSQDs
with narrow core band gap that encapsulated with a shell layer with much wider
band gap. To imply a strong confinement case, the band alignment is chosen
such that $V_{c}\gg E_{g}^{bulk}$, and $V_{v}\approx\frac{1}{2}E_{g}^{bulk}$.
For this purpose, we choose PbS/CdS core-shell QDs with material parameters as
shown in Table 1.

Table 1: List of parameters for PbS/CdS core-shell QDs. Materials | PbS [23] | CdS [23, 24]
---|---|---
Band gap (eV) | 0.41 | 2.5
$m_{e}/m_{0}$ | 0.08 | 0.2
$m_{h}/m_{0}$ | 0.08 | 0.7
Relative permittivity, $\varepsilon_{r}$ | 169 | 8.73
$V_{c}$ (eV) | 2.29
$V_{v}$ (eV) | 0.20

Using Eqs. (5), (8) and (9), the energy of electrons ($E_{e}$), holes
($E_{h}$) and the total transition energy in the type-II CSQDs
($E_{\text{qd}}$), we present their behaviors in Figure 2.

Figure 2: (color online) Energy of electrons (left panel) and holes (right
panel) for PbS/CdS. Transition energy for PbS/CdS CSQDs (bottom panel).

###  3.2 Type-II CSQDs with weak confinement

We apply our model single-band model to a wide core QDs, encapsulated with a
shell layer with slightly wider band gap. Then, to insinuate a weak
confinement case, the step-potential is chosen such that
$V_{c}\gtrapprox\frac{1}{2}E_{g}^{bulk}$, while keeping
$V_{v}\approx\frac{1}{2}E_{g}^{bulk}$. Here we choose ZnTe/ZnSe core-shell QDs
with material parameters as shown in Table 2.

Table 2: List of parameters for ZnTe/ZnSe core-shell QDs. Materials | ZnTe [5, 23] | ZnSe [5, 23]
---|---|---
Band gap (eV) | 2.394 | 2.8215
$m_{e}/m_{0}$ | 0.11 | 0.14
$m_{h}/m_{0}$ | 0.7 | 0.6
Relative permittivity, $\varepsilon_{r}$ | 8.7 | 9.1
$V_{c}$ (eV) | 1.6275
$V_{v}$ (eV) | 1.2000

As before from Eqs (5), (8) and (9), the energy of electrons ($E_{e}$), holes
($E_{h}$) and the total transition energy in the type-II CSQDs
($E_{\text{qd}}$) for weak confinement case are presented in Figure 3.

Figure 3: (color online) Energy of electrons (left panel) and holes (right
panel) for ZnTe/ZnSe. Transition energy for ZnTe/ZnSe CSQDs (bottom panel).

##  4 Radial probability

###  4.1 Strong confinement case

The overall trend for the transition energy in PbS/CdS CSQDs shows a jump
around 5$\%$ compared to bare PbS quantum dots when we added 0.25 nm of shell
thickness. However, it started to be decreasing after thicker shell layer
applied, as shown in bottom panel of Figure 2. This proves that capping a
quantum dots into a core-shell form will cause stronger localization while
inducing stress onto the charge carriers that will cause further lattice
mismatch [25]. The induced lattice mismatch now responsible to minimize the
strain caused by the strongly confined carriers. However, the trend of energy
of electrons keeps constant although the shell thickness was increased for any
given core radii. This is most probably due to the fact that the electrons is
now in the state of highly localized, which makes them to be less “triggered”
by the change in shell thickness. This can be proven by comparing the radial
probability of electrons in Figures 4 and 5 where different shell thickness
did not affect on how the electrons are localized.

Figure 4: (color online) Radial probability of carriers in 4.25 nm PbS/CdS
with 0.25 nm shell thickness. Figure 5: (color online) Radial probability of
carriers in 4.50 nm PbS/CdS with 0.5 nm shell thickness.

As for the holes, a jump in energy was shown in right panel of Figure 2,
before they started decreasing after 2.5 nm of CdS shell. As the strength of
confinement is much less in the valence band compared to in conduction band,
the holes are less likely to be highly localized as their electrons
counterpart, making them to be “aware” of any change in shell thickness.
Hence, as thickness is increased, the holes started to tunnel through between
core and shell regions, as shown by their spread of radial probability in
Figure 5. It can be seen that the spread of radial wavefunction of holes are
higher in Figure 5 compared to holes in 4.

###  4.2 Weak confinement case

The total transition energy in ZnTe/ZnSe CSQDs is higher compared to their
bare ZnTe quantum dots, by almost 2$\%$, as shown in bottom panel of Figure 3.
This lower jump is due to a weaker confinement strength compared to the case
of PbS/CdS in previous section. The trend of energy of electrons in ZnTe/ZnSe
is the same as in PbS/CdS, where their localized state causing them to be less
influenced by the increasing of shell thickness. Unlike in the case of strong
confinement where different core radii had different energy of holes, in weak
confinement, the difference in energy of holes are hardly noticeable. The
reason for this can be studied by observing their radial probability.

Although the ZnTe/ZnSe CSQDs has type-II band alignments, it mimics the
behaviour of type-I CSQDs where both charge carriers are likely to be confined
within the same region. It can be seen in Figures 6 and 7 that their holes are
likely (not entirely) to be localized within the core region instead of within
the shell band like any typical type-II CSQDs. From these results, it can be
inferred that weak confinement of the type-II case can be referred as pseudo-
type-II CSQDs [12].

Figure 6: (color online) Radial probability of carriers in 4.25 nm ZnTe/ZnSe
with 0.25 nm shell thickness. Figure 7: (color online) Radial probability of
carriers in 4.50 nm ZnTe/ZnSe with 0.5 nm shell thickness.

###  4.3 Effect of core size

The effects of different core sizes is further studied by analysing their
radial probabilities with the same shell thickness. In the case of strong
confinement of PbS/CdS CSQDs, we can see from bottom panel of Figure 2 that
the smallest CSQDs has the highest energy compared to the bigger ones. As
expected, this small size of CSQDs contributes to a greater localization of
charge carriers as shown in Figure 8, where carriers inside 4.00 nm core size
has the highest probability of being at the centre of the core.

Figure 8: (color onine) Radial probability of electrons (left panel) and holes
(right panel) in PbS/CdS CSQDs with 0.5 nm of shell thickness (strong
confinement case).

However, in the case of weak confinement of ZnTe/ZnSe CSQDs, the shell
thickness only influence the electrons confinement, and not the holes. We also
found that with different core radii, the holes still mostly spread throughout
the core region as observed in Figure 9. Hence, both Figure 7 and 9 confirm
that core and shell sizes will not influence on how the holes will be
localized or confined.

Figure 9: (color onine) Radial probability of electrons (left panel) and holes
(right panel) in ZnTe/ZnSe CSQDs with 0.5 nm of shell thickness weak
confinement case).

##  5 Conclusion

We have properly defined the strength of confinement based on their step-
potential rather than their size. The new conditions of strength confinement
was justified by investigating the effects of strong and weak confinement
cases in the type-II CSQDs with varying both core radius and shell thickness.
It can be seen that by introducing strong confinement in CQDS, the transition
energy would increase as high as 5$\%$, compared to the weak confinement case,
where the increase of their transition energy is only around 2$\%$. From the
theoretical point of view, this results demonstrate the flexibility to
manipulate the combination of core and shell materials while addressing the
strong and weak confinement conditions. Furthermore, the weakly confined
charge carriers in type-II CSQDs has the ability to demonstrate red-shift
emission where the transition energy of this CSQDs is lower compared to both
of ZnTe and ZnSe band gaps. This is favourable in fabricating near-infrared or
infrared photons emission with longer carriers lifetime [8].

Another worthy finding is how weak confinement case in type-II mimics
confinement behaviour of type-I CSQDs where both charge carriers are confined
within the same region, with their holes being less influenced by the changes
in both core and shell sizes. This has brought us to refer such case as
pseudo-type-II CSQDs. It is hoped that our findings can help experimentalist
to further understand the reason behind behavioural changes in type-II core-
shell quantum dots, which will help them to explain and justify their
experimental findings.

## Acknowledgments

We would like to thank Prof Zainal Abidin Talib for his kind guidance in
preparing this manuscript. This research was supported by Ministry of Higher
Education Malaysia (MOHE) through the Fundamental Research Grant Scheme
(FRGS/1/2018/STG02/UPM/02/10).

## References

  * [1] G. C. Dousmanis, Phys. Rev. Lett. 1, 55 (1958).
  * [2] M. Bruchez, M. Moronne, P. Gin, S. Weiss, and A. P. Alivisatos, Science 281, 2013 (1998); W. C. Chan and S. Nie, Science 281, 2016 (1998).
  * [3] X. Michalet, F.F. Pinaud, L.A. Bentolila, J.M. Tsay, S. Doose, J.J. Li, G. Sundaresan, A.M. Wu, S.S. Gambhir, and S. Weiss, Science 307, 538 (2005); X. Gao, Y. Cui, R.M. Levenson, L.W. Chung, and S. Nie, Nat. Biotechnol. 22, 969 (2004); W. Cai, D.W. Shin, K. Chen, O. Gheysens, Q. Cao, S.X. Wang, S.S. Gambhir, and X. Chen, Nano Lett. 6, 669 (2006).
  * [4] L. Shi, B. Hernandez, and M. Selke, J. Am. Chem. Soc. 128, 6278 (2006).
  * [5] G. Jia, B. Hao, X. Lu, and J. Yao, Int. J. Electrochem. Sci. 8, 8167 (2013).
  * [6] F. Poulsen and T. Hansen, J. Phys. Chem. C 121, 13655 (2017).
  * [7] N. Ca, N. Hien, N. Luyen, V. Lien, L. Thanh, P. Do, N. Bau, and T. Pham, J. Alloys Compd. 787, 823 (2019).
  * [8] S. Kim, B. Fisher, H.J. Eisler, and M. Bawendi, J. Am. Chem. Soc. 125, 11466 (2003).
  * [9] C. Vantakhah and A. Ebadi, Res. J. Recent Sci. 2, 21 (2013).
  * [10] R. Kostic and D. Stojanovic, J. Optoelectron. Adv. M. 6, 121 (2012).
  * [11] T. Cheche, V. Barna, and I. Stamatin, J. Optotoelectron. Adv. M. 15, 615 (2013).
  * [12] Y. Nandan and M. S. Mehata, Sci. Rep. 9, 2 (2019).
  * [13] D. Bimberg, M. Grundmann, and N. N. Ledentsov, Quantum dot heterostructures (John Wiley & Sons, 1999).
  * [14] D. Dorfs, T. Franzl, R. Osovsky, M. Brumer, E. Lifshitz, T. A. Klar, and A. Eychmuller, Small 4, 1148 (2008).
  * [15] B. Rocha, R. Jakomin, R. Kawabata, L. Dornelas, M. Pires, and P. Souza, 34th Symposium on Microelectronics Technology and Devices (SBMicro), IEEE (2019), pp. 1–3.
  * [16] Y. Shoji, R. Tamaki, and Y. Okada, AIP Advances 7, 065305 (2017).
  * [17] S. Kim, Y. T. Lim, E. G. Soltesz, A. M. De Grand, J. Lee, A. Nakayama, J. A. Parker, T. Mihaljevic, R. G. Laurence, D. M. Dor, et al., Nat Biotechnol 22, 93 (2004).
  * [18] T. Shelawati, M.S. Nurisya, A. Mazliana, et al., Superlattices and Microstructures 131, 95 (2019).
  * [19] O. Stier, M. Grundmann, and D. Bimberg, Phys. Rev B 59, 5688 (1999).
  * [20] A. L. Kaledin, D. Kong, K. Wu, T. Lian, and D. G. Musaev, J. Phys. Chem. C 122, 1874 (2018).
  * [21] E. Tyrrell and J. Smith, Phys. Rev. B 84, 165328 (2011).
  * [22] L. Brus, J. Phys. Chem. 90, 2555 (1986).
  * [23] U. K. Mishra and J. Singh, Semiconductor Device Physics and Design (Springer Netherlands, 2008).
  * [24] O. Madelung, U. Rossler, and M. Schulz, Semiconductors II-VI and I-VII Compounds; Semimagnetic Compounds: Cadmium sulfide (CdS) dielectric constants (Springer Berlin Heidelberg, Berlin, Heidelberg, 1999).
  * [25] K. Gong and D. F. Kelley, J. Chem. Phys. 141, 194704 (2014).