# Exact solution of a lambda quantum system driven by a two-photon wavepacket

W. Lopes 1 D. Valente 1<EMAIL_ADDRESS>1 Instituto de Física,
Universidade Federal de Mato Grosso, Cuiabá, MT, Brazil

###### Abstract

Three-level atoms in lambda configuration find diverse applications in quantum
information processing, and a promising way to manipulate their quantum states
is with single-photon pulses propagating in a waveguide (which can be
theoretically regarded as a highly broadband regime of the Jaynes-Cummings
model). Here, we analytically find the non-perturbative dynamics of a lambda
atom driven by a two-photon wavepacket, propagating in a one-dimensional
electromagnetic environment. As an application, we study the dynamics of a
quantum state purification. By comparing our exact model with an approximated
model of two cascaded single-photon wavepackets, we show how two-photon
nonlinearities and stimulated emission affect the purification.

##  I Introduction

Two-level quantum systems can play the role of qbits, i.e., quantum bits of
information. Implementing a logical qbit within the two ground states of a
three-level atom in $\Lambda$ configuration has the advantage that the excited
state provides an energetic barrier, blocking unwanted jumps between the two
target states, thus increasing the stability of physical operations. This
motivates various applications of the $\Lambda$ atoms, such as in quantum
memories rmp2017 ; molmer , quantum cloning dvnjp ; dv2012 ; qsr ; csimon1 ;
csimon2 , and quantum state purification qda ; cp22 .

To achieve efficient manipulation of the atom, one can either employ optical
cavities or one-dimensional (1D) waveguides. Waveguides are especially
desirable if one wishes that a propagating single-photon wavepacket interacts
with a single atom, for instance if the photon works as a flying qubit making
distant atoms in a network to communicate kimble ; lodahl ; nori ; jmg . From
a theoretical viewpoint, these engineered 1D environments can be described as
the highly broadband regime of the Jaynes-Cummings (JC) model chen2004 ; fan ;
domokos . This has unraveled a rich diversity of effects, including non-
Markovianity nonmarkov , broadband photon-photon correlations roy ; njp17 ;
lodahl22 , irreversible stimulated emission dvnjp ; continuum2018 ; dv2012 ,
as well as the application of quantum mirrors in Fabry-Perot prlfratini  and
Mach-Zehnder nelson  interferometers.

Here, our main motivation is to understand how the type of quantum state
purification for a $\Lambda$ atom as proposed in Ref.qda  is affected by a
two-photon wavepacket. To that end, we present an exact non-perturbative
solution, in real-space representation, for the dynamics of a two-photon
wavepacket interacting with a single $\Lambda$ system in a broadband one-
dimensional environment. We describe light-matter interaction by employing the
broadband JC model, and analytically compute the quantum state of the
composite atom-plus-field system within the two-excitation subspace. We obtain
the transition probability between the two ground states of the atom,
depending on the initial two-photon pulse parameters. In order to highlight
the presence of two-photon nonlinearities, we compare our exact results with
the effective model of two cascaded (consecutive) single-photon wavepackets.

Figure 1:  Quantum state purification with a two-photon wavepacket in a one-
dimensional environment with spatial coordinate $r$. (A) At time $t=0$, a two-
photon wavepacket $|2_{AA}\rangle$, with $AA$-polarized photons, both resonant
with transition frequency $\omega_{a}$, starts interacting with the atom. The
photons are described by their linewidths (inversely related with the spatial
spread of the pulse). The atom is initially in a generic mixture of
$|a\rangle$ and $|b\rangle$ states. The atom decay rates are $\Gamma_{a}$ and
$\Gamma_{b}$, respectively between states $|e\rangle$ to $|a\rangle$, and
$|e\rangle$ to $|b\rangle$. (B) Asymptotically, $t\rightarrow\infty$, the atom
may have been left at the pure state $|b\rangle$, as governed by the
transition probability $p_{a\rightarrow b}(\infty)$. The quantum state of the
outgoing wavepacket is altered to $|2_{AB}\rangle$, both in spacetime and in
polarization (in this case, to $AB$). Our broadband JC model, within a real-
space representation of the amplitudes, allows us to compute the time-
dependent probability amplitudes for the global atom-plus-field system.

##  II Model

We consider a global time-independent Hamiltonian of the system and its
environment,

$H=H_{S}+H_{I}+H_{E}.$ (1)

We model light-matter interaction via a dipole coupling, in the rotating-wave
approximation domokos ; chen2004 ; fan ; dvnjp ; dv2012 ; nonmarkov ,

$H_{I}=-i\hbar\sum_{\omega}(g_{a}a_{\omega}\sigma_{a}^{\dagger}+g_{b}b_{\omega}\sigma_{b}^{\dagger}-\mbox{H.c.}),$
(2)

where $\sigma_{k}=|k\rangle\langle e|$ (for $k=a,b$) and H.c. is the Hermitian
conjugate. Orthogonal polarization modes $\left\\{a_{\omega}\right\\}$ and
$\left\\{b_{\omega}\right\\}$ respectively interact with transitions
$|a\rangle$ to $|e\rangle$ and $|b\rangle$ to $|e\rangle$. This is a multi-
mode Jaynes-Cummings model, and our interest is in the limit of a broadband
continuum of frequencies,

$\sum_{\omega}\rightarrow\int d\omega\varrho_{\omega}\approx\varrho\int
d\omega,$ (3)

where we have defined a constant density of modes
$\varrho_{\omega}\approx\varrho$. Physically, this means that there is no
frequency filter in our idealized environment. The continuum of frequencies
allows us to employ a Wigner-Weisskopf approximation to obtain dissipation
rates $\Gamma_{a}=2\pi g_{a}^{2}\varrho$ and $\Gamma_{b}=2\pi
g_{b}^{2}\varrho$. The environment Hamiltonian reads

$H_{E}=\sum_{\omega}\hbar\omega(a_{\omega}^{\dagger}a_{\omega}+b_{\omega}^{\dagger}b_{\omega}).$
(4)

The system is described by a three-level Hamiltonian, with states named
$|a\rangle$, $|b\rangle$, and $|e\rangle$, so that

$H_{S}=\hbar\omega_{a}|e\rangle\langle e|+\hbar\delta_{ab}|b\rangle\langle
b|.$ (5)

When in $\Lambda$ configuration, these energy levels provide an elementary
energy barrier (the excited state $|e\rangle$) separating two highly stable
states (the two ground states $|a\rangle$ and $|b\rangle$). The energy
difference between states $|b\rangle$ and $|a\rangle$, given by
$\hbar\delta_{ab}=\hbar(\omega_{a}-\omega_{b})$, although relevant for thermal
equilibrium considerations, does not affect our results as far as
nonequilibrium effects are concerned.

As stated above, our main interest here is in applying our exact solution to
the study of quantum state purification. This is the process where the reduced
state of the atom asymptotically goes from a mixed to a pure state, such as

$p_{a}^{(0)}|a\rangle\langle a|+p_{b}^{(0)}|b\rangle\langle
b|\rightarrow|b\rangle\langle b|,$ (6)

after being driven by the photons and dissipated to the environment. To
achieve this kind of dynamics, we focus on the case where the initial state of
the field is a two-photon wavepacket polarized at modes
$\left\\{a_{\omega}\right\\}$. This means that the photons will leave the atom
unaltered, if initially at $|b\rangle$. Thus the reduced state of the atom
evolves as

$p_{a}^{(0)}\mbox{tr}_{F}[U|a,2_{AA}\rangle\langle
a,2_{AA}|U^{\dagger}]+p_{b}^{(0)}\mbox{tr}_{F}[U|b,2_{AA}\rangle\langle
b,2_{AA}|U^{\dagger}]=p_{a}^{(0)}\mbox{tr}_{F}[U|a,2_{AA}\rangle\langle
a,2_{AA}|U^{\dagger}]+p_{b}^{(0)}|b\rangle\langle b|,$ (7)

where $\mbox{tr}_{F}$ is the partial trace over the field degrees of freedom,
and $U=\exp(-iHt/\hbar)$. We thus have to focus on the
$|a\rangle\rightarrow|b\rangle$ transition probability, defined as

$p_{a\rightarrow b}(t)\equiv\langle b|\mbox{tr}_{F}[U|a,2_{AA}\rangle\langle
a,2_{AA}|U^{\dagger}]|b\rangle,$ (8)

if we wish to study the purification dynamics of the atomic state. If
$p_{a\rightarrow b}(\infty)=1$, the atom will be left at the pure state
$|b\rangle$, asymptotically. See Fig.(1).

##  III Results

###  III.1 Quantum state of the global atom-plus-field system in the two-
excitation subspace

We take advantage of the number conservation as provided by the JC model in
order to restrict ourselves to the two-excitation subspace, spanned by

$\displaystyle|\psi\rangle=\sum_{\omega}\psi_{\omega}^{A}\hat{a}_{\omega}^{\dagger}|e,0\rangle+\sum_{\omega}\psi_{\omega}^{B}\hat{b}_{\omega}^{\dagger}|e,0\rangle+\sum_{\omega_{1}\omega_{2}}\phi_{\omega_{1}\omega_{2}}^{AA}\hat{a}_{\omega_{1}}^{\dagger}\hat{a}_{\omega_{2}}^{\dagger}|a,0\rangle+2\sum_{\omega_{1}\omega_{2}}\phi_{\omega_{1}\omega_{2}}^{AB}\hat{a}_{\omega_{1}}^{\dagger}\hat{b}_{\omega_{2}}^{\dagger}|a,0\rangle$
$\displaystyle+\sum_{\omega_{1}\omega_{2}}\phi_{\omega_{1}\omega_{2}}^{BB}\hat{b}_{\omega_{1}}^{\dagger}\hat{b}_{\omega_{2}}^{\dagger}|b,0\rangle+2\sum_{\omega_{1}\omega_{2}}\phi_{\omega_{1}\omega_{2}}^{BA}\hat{b}_{\omega_{1}}^{\dagger}\hat{a}_{\omega_{2}}^{\dagger}|b,0\rangle.$
(9)

By substituting that into the Schrödinger equation,
$i\hbar\partial_{t}|\psi\rangle=H|\psi\rangle$, where $H$ is given in Eq.(1),
we find that

$\partial_{t}\psi_{\omega}^{A}(t)=-i(\omega_{a}+\omega)\psi_{\omega}^{A}(t)-2\sum_{\omega^{\prime}}[g_{a}\phi_{\omega\omega^{\prime}}^{AA}(t)+g_{b}\phi_{\omega\omega^{\prime}}^{BA}(t)].$
(10)
$\partial_{t}\psi_{\omega}^{B}(t)=-i(\omega_{a}+\omega)\psi_{\omega}^{B}(t)-2\sum_{\omega^{\prime}}[g_{b}\phi_{\omega\omega^{\prime}}^{BB}(t)+g_{a}\phi_{\omega\omega^{\prime}}^{AB}(t)];$
(11)
$\partial_{t}\phi_{\omega\omega^{\prime}}^{AA}(t)=-i(\omega+\omega^{\prime})\phi_{\omega\omega^{\prime}}^{AA}(t)+\frac{1}{2}(g_{a}\psi_{\omega}^{A}(t)+g_{a}\psi_{\omega^{\prime}}^{A}(t))$
(12)
$\partial_{t}\phi_{\omega\omega^{\prime}}^{BB}(t)=-i(\omega+\omega^{\prime}+\delta_{ab})\phi_{\omega\omega^{\prime}}^{BB}(t)+\frac{1}{2}(g_{b}\psi_{\omega}^{B}(t)+g_{b}\psi_{\omega^{\prime}}^{B}(t));$
(13)
$\partial_{t}\phi_{\omega\omega^{\prime}}^{AB}(t)=-i(\omega+\omega^{\prime})\phi_{\omega\omega^{\prime}}^{AB}(t)+\frac{1}{2}g_{a}\psi_{\omega^{\prime}}^{B}(t);$
(14)
$\partial_{t}\phi_{\omega\omega^{\prime}}^{BA}(t)=-i(\omega+\omega^{\prime}+\delta_{ab})\phi_{\omega\omega^{\prime}}^{BA}(t)+\frac{1}{2}g_{b}\psi_{\omega}^{A}(t).$
(15)

###  III.2 Formal solution of the amplitudes

Inspired by the Wigner-Weisskopf theory of the spontaneous emission, we first
integrate the field equations,

$\displaystyle\phi_{\omega\omega^{\prime}}^{BA}(t)=\phi_{\omega\omega^{\prime}}^{BA}(0)e^{-i(\omega+\omega^{\prime}+\delta_{ab})t}+\frac{1}{2}e^{-i(\omega+\omega^{\prime}+\delta_{ab})t}\int_{0}^{t}e^{i(\omega+\omega^{\prime}+\delta_{ab})t^{\prime}}g_{b}\psi_{\omega}^{A}(t^{\prime})dt^{\prime},$
(16)

and

$\displaystyle\phi_{\omega\omega^{\prime}}^{AA}(t)=\phi_{\omega\omega^{\prime}}^{AA}(0)e^{-i(\omega+\omega^{\prime})t}+\frac{1}{2}e^{-i(\omega+\omega^{\prime})t}\int_{0}^{t}e^{i(\omega+\omega^{\prime})t^{\prime}}\left(g_{a}\psi_{\omega}^{A}(t^{\prime})+g_{a}\psi_{\omega^{\prime}}^{A}(t^{\prime})\right)dt^{\prime}.$
(17)

We then substitute them back into the excited-state amplitudes,

$\displaystyle\frac{\partial}{\partial t}\psi_{\omega}^{A}(t)$
$\displaystyle=$
$\displaystyle-i(\omega_{a}+\omega)\psi_{\omega}^{A}(t)-2\sum_{\omega^{\prime}}[g_{a}\phi_{\omega\omega^{\prime}}^{AA}(0)e^{-i(\omega+\omega^{\prime})t}+g_{b}\phi_{\omega\omega^{\prime}}^{BA}(0)e^{-i(\omega+\omega^{\prime}+\delta_{ab})}+$
$\displaystyle-$ $\displaystyle
g_{a}\frac{1}{2}\int_{0}^{t}e^{-i(\omega+\omega^{\prime})(t-t^{\prime})}g_{a}\psi_{\omega}^{A}(t^{\prime})dt^{\prime}-g_{a}\frac{1}{2}\int_{0}^{t}e^{i(\omega+\omega^{\prime})(t-t^{\prime})}g_{a}\psi_{\omega^{\prime}}^{A}(t^{\prime})dt^{\prime}+$
$\displaystyle-$ $\displaystyle
g_{b}\frac{1}{2}\int_{0}^{t}e^{i(\omega+\omega^{\prime}+\delta_{ab})(t-t^{\prime})}g_{b}\psi_{\omega}^{A}(t^{\prime})dt^{\prime}].$
(19)

Because we are interested in the case where the initial state of the atom is
$|a\rangle$ and the field is prepared with two photons at polarization $A$, we
have only shown the solutions that are relevant to this scenario. The more
general case can be easily obtained by performing the same steps as used here.

###  III.3 Real-space representation

We define the real-space representations of all the amplitudes, namely,

$\displaystyle\psi^{A}(r,t)=\sum_{\omega}\psi_{\omega}^{A}(t)e^{ik_{\omega}r},$
(20)

and

$\displaystyle\phi^{AA}(r_{1},r_{2},t)=\sum_{\omega,\omega^{\prime}}\phi_{\omega,\omega^{\prime}}^{AA}(t)e^{ik_{\omega}r_{1}+ik_{\omega^{\prime}}r_{2}},$
(21)

with analogous expressions for the $B$ components. Here, we use a linear
approximation for the dispersion relation, $k_{\omega}=\omega/c$, with $c$
being the speed of light. The idea that only $k_{\omega}>0$ is considered can
be justified in two ways: either with unidirectional waveguides employing
chiral couplings between atom and field, or considering that each polarization
mode is a symmetric superposition of two counter-propagating modes in a
bidirectional waveguide.

In these representations, the equations of motion read

$\displaystyle\left[\frac{\partial}{\partial t}+c\frac{\partial}{\partial
r}\right]\psi^{A}(r,t)=-\left(\frac{\Gamma}{2}+i\omega_{a}\right)\psi^{A}(r,t)-2g_{a}\phi^{AA}(r-ct,-ct,0)$
$\displaystyle-\Gamma_{a}\psi^{A}(-r,t-r/c)\Theta(r)\Theta\left(t-\frac{r}{c}\right),$
(22)

where we have defined $\Gamma\equiv\Gamma_{a}+\Gamma_{b}$, as obtained after
we employ the Wigner-Weisskopf approximation,

$\displaystyle\Gamma_{a,b}=\sum_{\omega^{\prime}}2\pi
g_{a,b}^{2}\delta(\omega_{a,b}-\omega^{\prime})=2\pi g_{a,b}^{2}\varrho.$ (23)

Because the coupling $g_{a}$ is not experimentally accessible in this regime
(only $\Gamma_{a}$ is), we take $g_{a}=\sqrt{\Gamma_{a}/(2\pi\varrho)}$, where
$\varrho$ is the density of modes defined in Eq.(3) (the density of modes
$\varrho$ cancels out at the end, so no probabilities will depend on it).
Another crucial remark concerning the Wigner-Weisskopf approximation employed
here is that, at some point in the calculation, we have an integration of the
type $\int_{0}^{t}dt^{\prime}\delta(t^{\prime}-\tau)$, where $\tau=t-r/c$.
Because $0<t^{\prime}<t$, this term amounts to zero, unless $0<\tau<t$, which
gives us $0<t-r/c<t$; hence $0<t-r/c$ justifies the Heaviside step function
$\Theta(t-r/c)$, whereas $t-r/c<t$ gives us $r>0$, so $\Theta(r)$ appears. In
Eq.(22), we have also set $\phi_{\omega\omega^{\prime}}^{BA}(r-ct,-ct,0)=0$,
since we are considering that no $B$-polarized photons are incident on the
atom.

The interpretation of Eq.(22) is that of a wave equation propagating towards
the positive direction in space. The damping term accounts for the spontaneous
emission of the atom, so the homogeneous solution is
$\psi^{A}(r-ct,0)\exp[-(\Gamma+i\omega_{a})t]$. The second term on the right-
hand side of Eq.(22) corresponds to a pair of photons driving the atom, with
strength proportional to the coupling frequency $g_{a}$. If $g_{a}=g_{b}=0$,
the field freely propagates at speed $c$, so its amplitude will obey to
$\phi^{AA}(r_{1},r_{2},t)=\phi^{AA}(r_{1}-ct,r_{2}-ct,0)$, in that case. This
allows us to understand the following reasoning. The atom is at position
$r_{a}=0$, so if there is a chance that a photon starts at position
$r_{2}=-ct$, then this photon will reach the atom after a time $t$. This
explains the dependence of the driving amplitude on $-ct$. The other photon
started at $r_{1}=r-ct$, explaining why $\phi^{AA}(r-ct,-ct,0)$ also depends
on $r-ct$. Finally, the last term on the right-hand side describes stimulated
emission. To see that, we notice that $\Theta(r)$ is nonzero only for $r>0$,
so $\psi^{A}(-r,t-r/c)$ is the excited-state amplitude evaluated at a negative
position, $-r<0$, where the presence of a photon (with $A$ polarization) can
stimulate the emission of the excited atom.

###  III.4 Solution of the excited-state amplitude in the real-space
representation

Because Eq.(22) presents a discontinuity due to the Heaviside step functions,
we proceed to solve it by a stepwise partitioning the space into three
domains, as follows.

  * •

Region I: $r<0$.

$\displaystyle\left[\frac{\partial}{\partial t}+c\frac{\partial}{\partial
r}\right]\psi_{1}^{A}(r,t)=-\left(\frac{\Gamma}{2}+i\omega_{a}\right)\psi_{1}^{A}(r,t)-2g_{a}\phi^{AA}(r-ct,-ct,0).$
(24)

  * •

Region II: $0<r<ct$.

$\displaystyle\left[\frac{\partial}{\partial t}+c\frac{\partial}{\partial
r}\right]\psi_{2}^{A}(r,t)=-\left(\frac{\Gamma}{2}+i\omega_{a}\right)\psi_{2}^{A}(r,t)-\Gamma_{a}\psi_{1}^{A}(-r,t-r/c)$
$\displaystyle-2g_{a}\phi^{AA}(r-ct,-ct,0).$ (25)

Note that the dynamics in Region II depends on the solution of Region I.

  * •

Region III: $r>ct$.

$\displaystyle\left[\frac{\partial}{\partial t}+c\frac{\partial}{\partial
r}\right]\psi_{3}^{A}(r,t)=-\left(\frac{\Gamma}{2}+i\omega_{a}\right)\psi_{3}^{A}(r,t)-2g_{a}\phi^{AA}(r-ct,-ct,0).$
(26)

To solve for Region I, we look for a solution of the type

$\displaystyle\psi_{1}^{A}(r,t)=e^{-\left(\frac{\Gamma}{2}+i\omega_{a}\right)t}\bar{\psi}_{1}^{A}(r,t),$
(27)

as motivated by the homogeneous solution. We thus find that

$\displaystyle\hat{\theta}\bar{\psi}_{1}^{A}=-2g_{a}e^{\left(\frac{\Gamma}{2}+i\omega_{a}\right)t}\phi^{AA}(r-ct,-ct,0),$
(28)

where we have defined $\hat{\theta}\equiv\partial_{t}+c\partial_{r}$, so the
inhomogeneous solution is given by

$\displaystyle\bar{\psi}_{1}^{A}(r,t)=-2g_{a}\int_{0}^{t}e^{\left(\frac{\Gamma}{2}+i\omega_{a}\right)t^{\prime}}\phi^{AA}(r-ct,-ct^{\prime},0)dt^{\prime}.$
(29)

To make sense of this solution, we can see it as describing one photon driving
the atom (term $-ct^{\prime}$), while the other photon is located at position
$r<0$ at time $t>0$. A formal justification of the solution can be given by
assuming a general decomposition such as
$\phi^{AA}(r-ct,-ct,0)=\sum_{nm}C_{nm}R_{n}(r-ct)T_{m}(-ct)$, where $C_{nm}$
are generic coefficients, with $R_{n}$ and $T_{m}$ being generic functions. By
noting that $\hat{\theta}R_{n}(r-ct)=0$, we reobtain the original differential
equation.

Region II is far richer in terms of dynamical behavior, since it depends on
the solution of Region I. Again, we first eliminate the homogeneous part by
taking the transformation

$\psi_{2}^{A}(r,t)=e^{-\left(\frac{\Gamma}{2}+i\omega_{a}\right)t}\bar{\psi}_{2}^{A}(r,t),$
(30)

which implies that

$\displaystyle\hat{\theta}\bar{\psi}_{2}^{A}(r,t)=-e^{\left(\frac{\Gamma}{2}+i\omega_{a}\right)t}\left[\Gamma_{a}\psi_{1}^{A}\left(-r,t-r/c\right)+2g_{a}\phi^{AA}(r-ct,-ct,0)\right],$
(31)

where, from Eq.(29), we have that

$\displaystyle\psi_{1}^{A}(-r,t-r/c)$ $\displaystyle=$
$\displaystyle-2g_{a}e^{-\left(\frac{\Gamma}{2}+i\omega_{a}\right)(t-r/c)}\int_{0}^{t-r/c}e^{\left(\frac{\Gamma}{2}+i\omega_{a}\right)t^{\prime}}\phi^{AA}(-ct,-ct^{\prime},0)dt^{\prime}.$
(32) $\displaystyle\equiv$ $\displaystyle F(t,t-r/c).$ (33)

The purpose of defining $F(t,t-r/c)$ above is just to highlight that the
solution only depends on $t$ and $t-r/c$. Because both terms on the right-hand
side of Eq.(31) are functions of $t$ and $t-r/c$, we can search for a solution
in a similar form as we did in Region I, namely,

$\displaystyle\bar{\psi}_{2}^{A}(r,t)=$ $\displaystyle-$
$\displaystyle\Gamma_{a}\int_{t-r/c}^{t}e^{\left(\frac{\Gamma}{2}+i\omega_{a}\right)t^{\prime}}F(t^{\prime},t-r/c)dt^{\prime}-2g_{a}\int_{0}^{t}e^{\left(\frac{\Gamma}{2}+i\omega_{a}\right)t^{\prime}}\phi^{AA}(r-ct,-ct^{\prime},0)dt^{\prime}.$
(34)

A crucial remark: the integration in the first term of Eq.(34) is starting at
$t-r/c$. Replacing $t-r/c$ by $0$ in that term also provides a correct
solution for Eq.(31). However, our final interest is in solving Eq.(22), to
which this makes a key difference. We further explain why by the end of this
section.

Region III has precisely the same equation of motion of Region I, so the
solution is

$\displaystyle\bar{\psi}_{3}^{A}(r,t)=-2g_{a}e^{-\left(\frac{\Gamma}{2}+i\omega_{a}\right)t}\int_{0}^{t}e^{\left(\frac{\Gamma}{2}+i\omega_{a}\right)t^{\prime}}\phi^{AA}(r-ct,-ct^{\prime},0)dt^{\prime}.$
(35)

If one is interested in considering a two-photon wavepacket that initially
occupies only the $r<0$ region (Region I), as we are here, one can assume that
$\phi^{AA}(r_{1},r_{2},0)=0$ for $r_{1}>0$, which amounts to say that
$\phi^{AA}(r-ct,-ct^{\prime},0)=0$ for $r>ct$ (that defines Region III).

Putting all three regions together, and adding the homogeneous solution, we
get that

$\displaystyle\psi^{A}(r,t)$ $\displaystyle=$
$\displaystyle\psi^{A}(r-ct,0)e^{-\left(\frac{\Gamma}{2}+i\omega_{a}\right)t}-2g_{a}e^{-\left(\frac{\Gamma}{2}+i\omega_{a}\right)t}\int_{0}^{t}dt^{\prime}e^{\left(\frac{\Gamma}{2}+i\omega_{a}\right)t^{\prime}}\phi^{AA}(r-ct,-ct^{\prime},0)+$
(36) $\displaystyle+$ $\displaystyle
2\Gamma_{a}g_{a}\Theta(r)\Theta(t-r/c)e^{-\left(\frac{\Gamma}{2}+i\omega_{a}\right)t}e^{-\left(\frac{\Gamma}{2}+i\omega_{a}\right)(t-r/c)}$
$\displaystyle\times$
$\displaystyle\int_{t-r/c}^{t}dt^{\prime}e^{\left(\frac{\Gamma}{2}+i\omega_{a}\right)t^{\prime}}\int_{0}^{t-r/c}dt^{\prime\prime}e^{\left(\frac{\Gamma}{2}+i\omega_{a}\right)t^{\prime\prime}}\phi^{AA}(-ct^{\prime},-ct^{\prime\prime},0).$

The first term corresponds to two simultaneous events, namely, the free
evolution of the single-photon packet, and the spontaneous emission of the
atom. The second term corresponds to one photon exciting the atom, while the
other photon freely propagates. And the third term describes both photons
driving the atom, since $\phi^{AA}(-ct,-ct,0)$ indicates both photons having
arrived at the atomic position $r_{a}=0$ after time $t$ has passed.

We can now give a justification for the rather unusual integration of Eq.(34)
(from $t-r/c$ to $t$). Starting from Eq.(36) and substituting it back into
Eq.(22), we get a spurious extra term coming from
$\partial_{r}\Theta(r)=\delta(r)$. The lower limit defined as $t-r/c$ exactly
cancels out that spurious term, since $\delta(r)\int_{t-r/c}^{t}(...)=0$.

###  III.5 Solution of the field amplitudes in the real-space representation

Our main interest is in finding an expression for
$\phi^{BA}(r_{1},r_{2},t)=\sum_{\omega,\omega^{\prime}}\phi^{BA}_{\omega,\omega^{\prime}}(t)\exp(ik_{\omega}r_{1}+ik_{\omega^{\prime}}r_{2})$.
Analogous to Eq.(16), and choosing the initial condition
$\phi^{BA}_{\omega,\omega^{\prime}}(0)=0$, we have that

$\displaystyle\phi^{BA}(r_{1},r_{2},t)$ $\displaystyle=$
$\displaystyle\sum_{\omega\omega^{\prime}}\left(\frac{1}{2}e^{-i(\omega+\omega^{\prime}+\delta_{ab})t}\int_{0}^{t}e^{i(\omega+\omega^{\prime}+\delta_{ab})t^{\prime}}g_{b}\psi_{\omega}^{A}(t^{\prime})dt^{\prime}\right)e^{ik_{\omega}r_{1}+ik_{\omega^{\prime}}r_{2}}$
$\displaystyle=$
$\displaystyle\sum_{\omega\omega^{\prime}}\left(\frac{1}{2}\int_{0}^{t}e^{-i(\omega+\omega^{\prime}+\delta_{ab})(t-t^{\prime})}g_{b}\psi_{\omega}^{A}(t^{\prime})dt^{\prime}\right)e^{ik_{\omega}r_{1}+ik_{\omega^{\prime}}r_{2}}$
$\displaystyle=$
$\displaystyle\sum_{\omega\omega^{\prime}}\frac{1}{2}\int_{0}^{t}e^{ik_{\omega}[r_{1}-c(t-t^{\prime})]}e^{i\frac{\omega^{\prime}}{c}[r_{2}-c(t-t^{\prime})]}e^{-i\delta_{ab}(t-t^{\prime})}g_{b}\psi_{\omega}^{A}(t^{\prime})dt^{\prime}$
$\displaystyle=$
$\displaystyle\frac{1}{2}g_{b}\int_{0}^{t}\left(\sum_{\omega}\psi_{\omega}^{A}(t^{\prime})e^{ik_{\omega}[r_{1}-c(t-t^{\prime})]}\right)\sum_{\omega^{\prime}}e^{i\frac{\omega^{\prime}}{c}[r_{2}-c(t-t^{\prime})]}e^{-i\delta_{ab}(t-t^{\prime})}dt^{\prime}$
$\displaystyle=$
$\displaystyle\frac{1}{2}g_{b}\int_{0}^{t}\psi^{A}(r_{1}-c(t-t^{\prime}),t^{\prime})\sum_{\omega^{\prime}}e^{-i\omega^{\prime}[(t-r_{2}/c)-t^{\prime}]}e^{-i\delta_{ab}(t-t^{\prime})}dt^{\prime}.$

We take the continuum limit of the sum and perform a Wigner-Weisskopf
approximation once again. We also have to make sure that $0\leq t^{\prime}\leq
t$ is consistent with such approximation. These two inequalities are
guaranteed if we multiply the solution by $\Theta(r_{2})\Theta(t-r_{2}/c)$. We
finally find that

$\displaystyle\phi^{BA}(r_{1},r_{2},t)=\frac{\sqrt{2\pi\varrho\Gamma_{b}}}{2}\Theta(r_{2})\Theta\left(t-\frac{r_{2}}{c}\right)e^{-i\delta_{ab}r_{2}/{c}}\psi^{A}(r_{1}-r_{2},t-r_{2}/c).$
(37)

###  III.6 Transition probabilities in the real-space representation

The transition probability between the two ground states of the $\Lambda$
atom, i.e., from $|a\rangle$ to $|b\rangle$ is given by

$p_{a\rightarrow b}(t)=\langle
b|\mbox{tr}_{F}(|\psi(t)\rangle\langle\psi(t)|)|b\rangle,$ (38)

provided that
$|\psi(0)\rangle=|a,2_{AA}\rangle=\sum_{\omega_{1}\omega_{2}}\phi_{\omega_{1}\omega_{2}}^{AA}(0)\hat{a}_{\omega_{1}}^{\dagger}\hat{a}_{\omega_{2}}^{\dagger}|a,0\rangle$.
Also, $|\psi(t)\rangle=\exp(-iHt/\hbar)|\psi(0)\rangle$, where $H$ is given in
Eq.(1), and $tr_{F}$ is the partial trace over the field degrees of freedom.

Within the two-excitation subspace, Eq.(9), we have that

$p_{a\rightarrow
b}(t)=\sum_{\omega_{1}\omega_{2}}|\phi_{\omega_{1}\omega_{2}}^{BB}(t)|^{2}+4\sum_{\omega_{1}\omega_{2}}|\phi_{\omega_{1}\omega_{2}}^{BA}(t)|^{2}.$
(39)

Here, we are mainly interested in the case where
$\phi^{BB}_{\omega_{1},\omega_{2}}(0)=0$, implying that
$\phi^{BB}_{\omega_{1},\omega_{2}}(t)=0$. To compute the non-vanishing term,
we have to transform it to the real-space representation, finding that

$\sum_{\omega_{1}\omega_{2}}|\phi_{\omega_{1}\omega_{2}}^{BA}(t)|^{2}=\frac{1}{(2\pi\varrho
c)^{2}}\int_{-\infty}^{\infty}\int_{-\infty}^{\infty}|\phi^{BA}(r_{1},r_{2},t)\big{|}^{2}dr_{1}dr_{2}.$
(40)

By substituting Eq.(37) in the double integration above, and using Eq.(36) for
the excited-state amplitude, we have the exact solution of the transition
probability for a generic initial pulse $\phi^{AA}(r_{1},r_{2},0)$.

###  III.7 Analytically solvable example

Let us compute

$\displaystyle p_{a\rightarrow b}(t)$ $\displaystyle=\frac{4}{(2\pi\varrho
c)^{2}}\int_{-\infty}^{\infty}dr_{2}\int_{-\infty}^{\infty}dr_{1}\
|\phi^{BA}(r_{1},r_{2},t)\big{|}^{2}$ (41)
$\displaystyle=\frac{2\pi\varrho\Gamma_{b}}{(2\pi\varrho
c)^{2}}\int_{0}^{ct}dr_{2}\int_{-\infty}^{\infty}dr_{1}\
|\psi^{A}(r_{1}-r_{2},t-r_{2}/c)\big{|}^{2},$ (42)

for the case where the atom is initially at ground-state $|a\rangle$, so that
$\psi^{A}(r,0)=0$, and the initial photon pair is prepared by means of the
spontaneous emission of two hypothetical uncorrelated atoms, arbitrarily far
away from our $\Lambda$ atom, both polarized at $A$. We can thus choose

$\phi^{AA}(r_{1},r_{2},0)=N\bigg{[}f_{1}(r_{1})f_{2}(r_{2})+f_{1}(r_{2})f_{2}(r_{1})\bigg{]},$
(43)

where

$f_{k}(r)=\Theta(-r)\exp[(\Delta_{k}/2+i\omega_{Lk})r/c],$ (44)

and $N$ is a normalization constant satisfying
$\langle\psi(0)|\psi(0)\rangle=1$, equivalent to
$2\sum_{\omega_{1},\omega_{2}}|\phi^{AA}_{\omega_{1},\omega_{2}}(0)|^{2}=1$,
or in the real-space representation,

$\frac{2}{(2\pi\varrho
c)^{2}}\int_{-\infty}^{\infty}\int_{-\infty}^{\infty}|\phi^{AA}(r_{1},r_{2},0)\big{|}^{2}dr_{1}dr_{2}=1.$
(45)

Using Eqs.(43) and (44), and also assuming resonant photons,
$\omega_{L1}=\omega_{L2}$, we find that

$N=\pi\varrho\sqrt{\frac{\Delta_{1}\Delta_{2}}{1+r_{12}}},$ (46)

where we have defined the ratio $r_{12}\equiv
4\Delta_{1}\Delta_{2}/(\Delta_{1}+\Delta_{2})^{2}$. Note that $r_{12}\leq 1$,
the equality holding if $\Delta_{1}=\Delta_{2}$.

We now have two terms to solve in Eq.(36) for $\psi^{A}(r,t)$.

The first, where the integration of $\phi^{AA}(r-ct,-ct^{\prime},0)$ appears,
amounts to

$\psi^{A1}(r,t)=-2Ng_{a}\Theta(t-r/c)e^{-\frac{\Gamma}{2}t}\bigg{[}e^{-\frac{\Delta_{1}}{2}(t-\frac{r}{c})}\frac{2}{\Gamma-\Delta_{2}}\left(e^{\left(\frac{\Gamma-\Delta_{2}}{2}\right)t}-1\right)+e^{-\frac{\Delta_{2}}{2}(t-\frac{r}{c})}\frac{2}{\Gamma-\Delta_{1}}\left(e^{\left(\frac{\Gamma-\Delta_{1}}{2}\right)t}-1\right)\bigg{]}.$
(47)

When we compute $\psi^{A1}(r_{1}-r_{2},t-r_{2}/c)$, we get that
$\Theta(t-r_{1}/c)$ alters the integration from
$\int_{0}^{ct}dr_{2}\int_{-\infty}^{\infty}dr_{1}$ to
$\int_{0}^{ct}dr_{2}\int_{-\infty}^{ct}dr_{1}$.

The second term, where a double integration of
$\phi^{AA}(-ct^{\prime},-ct^{\prime\prime},0)$ appears, amounts to

$\displaystyle\psi^{A2}(r,t)$ $\displaystyle=$
$\displaystyle\frac{8N}{(\Gamma-\Delta_{1})(\Gamma-\Delta_{2})}g_{a}\Gamma_{a}\Theta(r)\Theta(t-r/c)e^{-\frac{\Gamma}{2}t}e^{-\frac{\Gamma}{2}\left(t-\frac{r}{c}\right)}$
(48) $\displaystyle\times$
$\displaystyle\bigg{[}\left(e^{\left(\frac{\Gamma-\Delta_{1}}{2}\right)t}-e^{\left(\frac{\Gamma-\Delta_{1}}{2}\right)\left(t-\frac{r}{c}\right)}\right)\left(e^{\left(\frac{\Gamma-\Delta_{2}}{2}\right)\left(t-\frac{r}{c}\right)}-1\right)+$
$\displaystyle+$
$\displaystyle\left(e^{\left(\frac{\Gamma-\Delta_{2}}{2}\right)t}-e^{\left(\frac{\Gamma-\Delta_{2}}{2}\right)\left(t-\frac{r}{c}\right)}\right)\left(e^{\left(\frac{\Gamma-\Delta_{1}}{2}\right)\left(t-\frac{r}{c}\right)}-1\right)\bigg{]}.$

In the terms proportional to $\psi^{A2}(r_{1}-r_{2},t-r_{2}/c)$, we have step
functions in the form $\Theta(r_{1}-r_{2})\Theta(t-r_{1}/c)$, altering the
integration from $\int_{0}^{ct}dr_{2}\int_{-\infty}^{\infty}dr_{1}$ to
$\int_{0}^{ct}dr_{2}\int_{r_{2}}^{ct}dr_{1}$.

We now have to compute, with the help of Eqs.(47) and (48),

$\displaystyle p_{a\rightarrow b}(t)$
$\displaystyle=\frac{2\pi\varrho\Gamma_{b}}{(2\pi\varrho
c)^{2}}\int_{0}^{ct}dr_{2}\int_{-\infty}^{ct}dr_{1}\big{|}\psi^{A1}(r_{1}-r_{2},t-r_{2}/c)\big{|}^{2}$
$\displaystyle+\frac{2\pi\varrho\Gamma_{b}}{(2\pi\varrho
c)^{2}}\int_{0}^{ct}dr_{2}\int_{r_{2}}^{ct}dr_{1}\big{|}\psi^{A2}(r_{1}-r_{2},t-r_{2}/c)\big{|}^{2}$
$\displaystyle+\frac{2\pi\varrho\Gamma_{b}}{(2\pi\varrho
c)^{2}}\int_{0}^{ct}dr_{2}\int_{r_{2}}^{ct}dr_{1}2\psi^{A1}(r_{1}-r_{2},t-r_{2}/c)\psi^{A2}(r_{1}-r_{2},t-r_{2}/c),$
(49)

where the Heaviside step functions have been already considered in the
integration limits.

Because the analytical expressions are too cumbersome, in Fig.(2) we plot
$p_{a\rightarrow b}(\infty)$ as a function of $\Delta_{2}$ (in log-linear
scale), for various $\Delta_{1}$ (blue curves). We consider both
$\Gamma_{b}=\Gamma_{a}$ and $\Gamma_{b}\neq\Gamma_{a}$ scenarios.

For the sake of a better understanding, we compare our exact two-photon
wavepacket model with the approximated model where two cascaded (consecutive,
uncorrelated, distinguishable) photons drive the atom (orange curves). In this
cascaded model, we have that the probability transition is given by ganascini

$p^{(C)}_{a\rightarrow b}(\infty)=p_{1}+(1-p_{1})p_{2},$ (50)

with $p_{k}=r_{\Gamma}/(1+\Delta_{k}/\Gamma)$, and $r_{\Gamma}\equiv
4\Gamma_{a}\Gamma_{b}/\Gamma^{2}$. Eq.(50) has a simple interpretation: either
the first photon drives the transition, with probability $p_{1}$, or it leaves
the atom unaltered, with probability $1-p_{1}$, and the second photon drives
the transition, with probability $p_{2}$. Note that $p^{(C)}_{a\rightarrow
b}(\infty)=1$ can only be obtained if $r_{\Gamma}=1$, which means that equal
decay rates, $\Gamma_{a}=\Gamma_{b}$, are a necessary condition for a full
purification to be achieved with two photons.

Figure 2:  Asymptotic transition probabilities as functions of the linewidth
$\Delta_{2}$ (in $\Gamma_{a}$ units). Exact model: $p_{a\rightarrow
b}(\infty)$ (blue). Approximated (cascaded) model: $p^{(C)}_{a\rightarrow
b}(\infty)$ (orange). We set $\Gamma_{a}=1$ in all plots. (A)
$\Gamma_{b}=0.5$, and $\Delta_{1}=0.001$. (B) $\Gamma_{b}=0.5$, and
$\Delta_{1}=100$. (C) $\Gamma_{b}=0.5$, and $\Delta_{1}=0.5$. (D)
$\Gamma_{b}=\Gamma_{a}=1$, and $\Delta_{1}=0.1$. Inset:
$\Gamma_{b}=\Gamma_{a}=1$, and $\Delta_{1}=0.001$. (A) and (B) show the
equivalences between the models, in the very low excitation regimes. (C) and
(D) show the deviations between the models, due to two-photon nonlinearities.
The dips (non-monotonic behaviors) in the blue curves of panels (C) and (D)
indicate the presence of stimulated emission.

Fig.(2)(A) and (B) show regimes where the exact ($p_{a\rightarrow b}(\infty)$,
blue) and the approximate ($p^{(C)}_{a\rightarrow b}(\infty)$, orange) models
are equivalent (blue and orange curves coincide). This happens either for too
small, $\Delta_{1}/\Gamma_{a}=0.001$ at (A), or too large linewidths,
$\Delta_{1}/\Gamma_{a}=100$ at (B). In both cases, the photons have very small
chance of exciting the atom, so two-photon nonlinearities are negligible.

Fig.(2)(C) and (D), by contrast, show how two-photon nonlinearities can
violate the cascaded-photons approximation for $\Delta_{2}/\Gamma_{a}\sim 1$.
In (C), we have set $\Gamma_{b}=\Gamma_{a}/2$, and $\Delta_{1}=\Gamma_{a}/2$,
thus finding a very pronounced deviation. In (D), we have set
$\Gamma_{b}=\Gamma_{a}$, and $\Delta_{1}/\Gamma_{a}=0.1$, showing a less
pronounced deviation, due to the fact that equal decay rates allow for much
closer-to-one transition probabilities then the unequal decay rates scenario.
To give further evidence of this, we also plot the case where
$\Gamma_{a}=\Gamma_{b}$, and $\Delta_{1}=0.001$, to show that both curves
saturate the transition probabilities, $p_{a\rightarrow b}(\infty)\approx
p^{(C)}_{a\rightarrow b}(\infty)\rightarrow 1$ for all values of $\Delta_{2}$
(inset).

It is important to emphasize the qualitative difference between the two
models, as shown in panels (C) and (D): the blue curve is non-monotonic,
unlike the orange curve. This behavior means that overlapping photons of
intermediate linewidths ($\Delta_{1}\sim\Delta_{2}\sim\Gamma$) can suppress
(rather than simply not contributing to) the transition from $|a\rangle$ to
$|b\rangle$. In other words, very large linewidths are unable to drive the
transition, but intermediate linewidths actively tend to preclude it to take
place. This negative two-photon nonlinearity effect can be understood in terms
of a stimulated emission (see discussions in Ref.dvnjp , for instance). If the
first photon has appreciable chance of exciting the atom, the presence of
another photon of similar linewidth can induce stimulated emission to take
place, forcing the atom to go back to initial state $|a\rangle$, instead of
transitioning to $|b\rangle$. This creates the dip in the blue curves of
panels (C) and (D).

##  IV Conclusions

We have developed a real-space representation of a two-photon wavepacket
driving a single three-level atom in $\Lambda$ configuration. We have solved
the non-perturbative dynamics, as modeled by the broadband JC model. We have
found analytical expressions both for the excited-state amplitudes (see
Eq.(36)) and the ground state amplitudes (see Eq.(37)). Because we were
interested in calculating the probability $p_{a\rightarrow b}(t)$ of
transition from $|a\rangle$ to $|b\rangle$, we have focused on $\phi^{BA}$, as
well as on $\psi^{A}$. Nevertheless, our method is directly applicable to
compute all the other amplitudes as defined within the two-excitation
subspace. We have also worked out an analytically solvable example, consisting
of an initial wavepacket of two photons of exponential envelope shape,
resonant with the atom $|a\rangle$ to $|e\rangle$ transition, and analyzed how
the probability transition behaves for diverse linewidths of the photons. We
have found that two-photon nonlinearities (taking place when both photons have
linewidths comparable to the atomic decay rates) can slightly suppress the
transition, due to the non-negligible effect of stimulated emission.

In future studies, we intend to explore further applications of our general
solutions. For instance, we would like to find how a two-photon wavepacket
affects quantum cloning processes, such as that from Ref.qsr . Also, we plan
to better understand the thermodynamics of this processes. How to define work
in this case, and how does the work relate with the probability transition
from $|a\rangle$ to $|b\rangle$? These would mean a generalization of the
quantum dissipative adaptation as proposed in Ref.qda .

###### Acknowledgements.

This work was supported by CNPq (402074/2023-8) and INCT-IQ (465469/2014-0),
Brazil. W.L. was supported by CAPES.

## References

  * (1) N. V. Vitanov, A. A. Rangelov, B. W. Shore, K. Bergmann, Stimulated Raman adiabatic passage in physics, chemistry, and beyond Rev. Mod. Phys. 89, 015006 (2017).
  * (2) K. Debnath, A. H. Kiilerich, A. Benseny, and K. Molmer, Coherent spectral hole burning and qubit isolation by stimulated Raman adiabatic passage, Phys. Rev. A 100, 023813 (2019).
  * (3) C. Simon, G. Weihs, and A. Zeilinger, Optimal quantum cloning via stimulated emission. Phys. Rev. Lett. 84, 2993 (2000).
  * (4) J. Kempe, C. Simon, and G. Weihs, Optimal photon cloning. Phys. Rev. A 62, 032302 (2000).
  * (5) D. Valente et. al., Optimal irreversible stimulated emission, New J. Phys. 14, 083029 (2012).
  * (6) Valente, D. et al. Universal optimal broadband photon cloning and entanglement creation in one-dimensional atoms. Phys. Rev. A 86, 022333 (2012).
  * (7) D. Valente, Self-replication of a quantum artificial organism driven by single-photon pulses, Sci. Rep. 11, 16433 (2021)
  * (8) D. Valente, F. Brito, and T. Werlang, Quantum dissipative adaptation, Comm. Phys. 4, 11 (2021).
  * (9) T. Werlang, M. Matos, F. Brito, and D. Valente, Emergence of energy-avoiding and energy-seeking behaviors in nonequilibrium dissipative quantum systems, Comm. Phys. 5, 7 (2022).
  * (10) H. J. Kimble, The quantum internet. Nature 453, 1023 (2008).
  * (11) J. Claudon et. al., A highly efficient single-photon source based on a quantum dot in a photonic nanowire. Nat. Photon. 4, 174–177 (2010).
  * (12) Lodahl, P., Mahmoodian, S., and Stobbe. S. Interfacing single photons and single quantum dots with photonic nanostructures. Rev. Mod. Phys. 87, 347 (2015).
  * (13) X. Gu., A. F. Kockum, A. Miranowicz _et al._ , Microwave photonics with superconducting quantum circuits. Phys. Rep. 718, 1-102 (2017).
  * (14) P. Domokos, P. Horak, and H. Ritsch, Quantum description of light-pulse scattering on a single atom in waveguides, Phys. Rev. A 65, 033832 (2002).
  * (15) Chen, T. W., Law, C. K., and Leung, P. T. Single-photon scattering and quantum-state transformations in cavity QED. Phys. Rev. A 69, 063810 (2004).
  * (16) J. T. Shen and S. Fan, Coherent photon transport from spontaneous emission in one-dimensional waveguides, Opt. Express 30, 2001 (2005).
  * (17) Valente, D., Arruda, M. F. Z., and Werlang, T. Non-Markovianity induced by a single-photon wave packet in a one-dimensional waveguide. Opt. Lett. 41, 3126-3129 (2016).
  * (18) Dibyendu Roy, C. M. Wilson, and O. Firstenberg Colloquium: strongly interacting photons in one-dimensional continuum. Rev. Mod. Phys. 89, 021001 (2017).
  * (19) Anders Nysteen et. al., Scattering of two photons on a quantum emitter in a one-dimensional waveguide: exact dynamics and induced correlations. New J. Phys. 17, 023030 (2015).
  * (20) A H. L. Jeannic, et. al., Dynamical photon–photon interaction mediated by a quantum emitter. Nat Phys 18, 1191–1195 (2022).
  * (21) K. A. Fischer, Exact calculation of stimulated emission driven by pulsed light. OSA Continuum 1, 772 (2018).
  * (22) F. Fratini et. al., Fabry-Perot Interferometer with Quantum Mirrors: Nonlinear Light Transport and Rectification. Phys. Rev. Lett. 113, 243601 (2014).
  * (23) N. Almeida, T. Werlang, and D. Valente, Mach-Zehnder interferometer with quantum beamsplitters, J. Opt. Soc. Am. B 36, 3357-3363 (2019).
  * (24) T. Ganascini, T. Werlang, and D. Valente, Quantum dissipative adaptation with cascaded photons. arXiv:2312.01157.


# Dealing with Cross-Task Class Discrimination in Online Continual Learning

Yiduo Guo1, Bing Liu4, Dongyan Zhao1,2,3
1Wangxuan Institute of Computer Technology, Peking University. 2BIGAI,
Beijing, China. 3National Key Laboratory of
General Artificial Intelligence. 4Department of Computer Science, University
of Illinois Chicago
<EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS>
###### Abstract

Existing continual learning (CL) research regards catastrophic forgetting (CF)
as almost the only challenge. This paper argues for another challenge in
class-incremental learning (CIL), which we call cross-task class
discrimination (CTCD), i.e., how to establish decision boundaries between the
classes of the new task and old tasks with no (or limited) access to the old
task data. CTCD is implicitly and partially dealt with by replay-based
methods. A replay method saves a small amount of data (replay data) from
previous tasks. When a batch of current task data arrives, the system jointly
trains the new data and some sampled replay data. The replay data enables the
system to partially learn the decision boundaries between the new classes and
the old classes as the amount of the saved data is small. However, this paper
argues that the replay approach also has a dynamic training bias issue which
reduces the effectiveness of the replay data in solving the CTCD problem. A
novel optimization objective with a gradient-based adaptive method is proposed
to dynamically deal with the problem in the online CL process. Experimental
results show that the new method achieves much better results in online CL.

##  1 Introduction

Continual learning (CL) learns a sequence of tasks incrementally. This work
focuses on the class incremental learning (CIL) setting [32] in online CL. In
CIL, each task consists of a set of unique classes, the sets of classes of any
two different tasks are disjoint and the system has no access to the task
information in testing. In online CL, the data comes gradually from a data
stream. Whenever the small batch of data arrives, it is trained in one
iteration. Thus, the data for each task is effectively trained in one epoch.

Existing CL papers almost regard catastrophic forgetting (CF) as the only
issue for CL. In fact, CIL also has another major challenge. When the system
learns a new task, if no data from previous tasks is available, it has no way
to establish decision boundaries between new classes and old classes in
previous tasks. Even if there is no CF, the classification results will still
be poor. We call this problem, cross-task class discrimination (CTCD). Those
approaches that do not save any previous data, e.g., regularization-based or
orthogonal projection-based, do not deal with CTCD. Replay-based methods
implicitly deal with CTCD to some extent because such a method uses a memory
buffer $\mathcal{M}$ to save a small amount of data (replay data) from old
tasks. When a small batch of current task data $X^{\textit{new}}$ arrives, the
system jointly trains $X^{\textit{new}}$ and some sampled replay data
$X^{\textit{buf}}$ from $\mathcal{M}$. $X^{\textit{buf}}$ enables the system
to partially learn the decision boundaries between the new classes and the old
classes because the amount of the saved data is very small.

Due to the limited replay data, the training is biased, which reduces its
ability to solve the CTCD problem. To make matters worse, the training bias
also changes as more tasks are learned. This paper first shows that the
problem is reflected as gradient imbalance (GI) on logits, i.e., higher
positive gradients than negative gradients on the logits and vice versa. It
further shows that GI is caused by two main issues. The first is data
imbalance. Since the memory buffer size, the batch size of the new data
$X^{\textit{new}}$, and the sampled data $X^{\textit{buf}}$ from the memory
buffer are all fixed, if the system has learned many tasks, the average number
of samples in each previous class in $X^{\textit{buf}}$ will be much smaller
than that of each class in $X^{\textit{new}}$. This results in higher positive
gradients than negative gradients on the logits of the previous classes
leading to training bias and poor decision boundaries (or weak CTCD
capability) between the classes of the new and old tasks. The second is CL
imbalance, i.e., CL training focuses more on the new samples (which are harder
to train as they are new) than the replayed samples (which have been seen and
trained many times before). This causes further GI. This imbalance is involved
(see Sec. 4.2 for details).

Some existing works[2, 42] have tried to deal with data imbalance in offline
CL. For example, SSIL[2] separately calculates the cross-entropy loss of the
new data and the replay data to mitigate data imbalance. But they are not from
the gradient angle. The second issue of GI is more complex and has not been
attempted before.

This paper proposes a novel method, called GSA (Gradient Self-Adaptation), to
deal with GI (and CTCD) in online CL. GSA includes a new training objective
and a gradient-based self-adaptive loss to compensate for the GI. The loss is
dynamically controlled by two gradient rates which automatically measure and
adapt to the dynamic GI situation. The main contributions of this paper are:

(1) It deals with the CTCD problem in online CL and proposes a new
optimization framework that decomposes the problem into cross-task
classification and within-task classification (see Section 5). In [22], CTCD
is called inter-task class separation, but it uses an out-of-distribution
based approach to dealing with the problem in offline CL. The paper uses a
replay-based approach for online CL.

(2) It analyzes the CTCD problem from the gradient imbalance (GI) angle and
finds two kinds of gradient imbalance (data imbalance and CL imbalance) (see
Section 4). Based on the analysis, it proposes a gradient-based self-adaptive
loss to compensate for the GI.

(3) Experiments in both the disjoint and long-tail online CL settings show
that GSA outperforms strong baselines by a large margin (see Section 6).

##  2 Related Work

There are many existing CL techniques. Regularization-based approaches
penalize changes to important parameters of old tasks [23, 50, 52]. Replay-
based approaches save some past data and replay them in new task training [12,
8, 37, 13, 54, 24, 47, 44]. Generative replay builds data generators to
generate pseudo old data for replaying [41, 38, 45, 21, 18, 34]. Dynamic
architectures based approaches [43, 39, 1, 49, 28, 20, 36, 14] overcomes CF by
expanding or isolating parameters. Data augmentation has also been used to
learn better features for CL recently [55].

Online CL and replay methods: Existing online CL methods mainly use the replay
approach. ER randomly samples the buffer data [11]. MIR chooses buffer data
whose loss increases the most [3]. Shapley value theory is applied in ASER for
memory update/retrieval [40]. Knowledge distillation is employed in DER++ [6].
NCCL deals with CF by calibrating the network [48]. Contrastive learning is
used in SCR [30]. GDumb samples and stores class-balanced data in memory [35].
However, it does not deal with CL imbalance identified in this paper. GSS
diversifies the gradients of the samples in the buffer [4]. ER-AML [7] avoids
CF by shielding the learned representations from drastic adaptation to
accommodate new classes. OCM is based on mutual information [16]. Our work
identifies two gradient imbalance issues and designed a new training strategy
and a novel loss function to deal with them.

Data imbalance in CL: Several researchers have dealt with data imbalance in
offline CL. SSIL [2] isolates the computation of softmax on previous and new
classes to update the model. LUCIR [19] uses the cosine normalization to
calculate the predicted probability and a margin ranking loss to separate
ground-truth old classes from new classes. [8] adds an additional fine-tuning
stage with a small learning rate and a balanced subset of samples. BiC [46]
adds a bias correction layer. [53] uses weight alignment to correct the biased
weights. CCIL [33] applies a loss to balance intra-task and inter-task
learning. However, many of these algorithms [53, 8] need the full data of the
current task to be available upfront and/or multiple training epochs to
address the data imbalance issue, which are not suitable for online CL as
online CL does not have the full training data of a task available when the
task arrives in the data stream.

Gradient-based CL methods: GEM [29] and A-GEM [10] rotate the current gradient
when the angle between the current gradient and the gradient computed on the
reference memory is obtuse. MEGA [17] uses the loss of $X^{\textit{new}}$ and
$X^{\textit{buf}}$ to adjust the relative importance of learning a new task
and maintaining the past knowledge. OWM [51], OGD [14], and AOP [15] avoid CF
by projecting the gradients on the new task onto an orthogonal subspace of old
tasks. These methods do not deal with gradient imbalance but represent a
different family of CL approaches.

##  3 Preliminary

Problem description. Following [3], we learn a sequence of tasks, 1, 2, …,
$t$, …, in the online CIL setting. In online CL, each incrementally arriving
data batch is seen only once by the system. So the system effectively learns
each task in one epoch. For a task $t$, we denote its dataset as
$\\{(x^{(t)}_{k},y^{(t)}_{k})\\}^{n^{t}}_{k=1}$, where $n^{t}$ is the number
of training samples in task $t$. $L^{(t)}$ denotes the set of classes of task
$t$.

Model architecture and training setting. Our model $F$ consists of a feature
extractor $h_{\theta}$ with the parameter set $\theta$ and a classifier
$f_{\phi}$ with the parameter set $\phi$. It uses a replay method with a
memory buffer $\mathcal{M}$. Whenever a small batch of new data
$X^{\textit{new}}$ from the data stream is accumulated, it is trained jointly
with a small batch of data $X^{\textit{buf}}$ sampled from $\mathcal{M}$ to
update the model in one training iteration. The model produces the logits
$F(x;\theta,\phi)=f_{\phi}(h_{\theta}(x))$, which are used to calculate the
loss or to predict in testing. Reservoir sampling is used for memory update.

Gradient imbalance (GI) on logits. We now introduce the proposed gradient-
based analysis on logits and the concept of GI. In CL, each task is usually
learned by minimizing the softmax cross-entropy loss, $\mathcal{L}_{ce}$:

$\displaystyle\mathcal{L}_{ce}(o(x;\theta,\phi))$
$\displaystyle=-\sum^{|C_{\textit{seen}}|}_{i=1}l_{c_{i}}\log(p^{c_{i}}),\qquad
p^{c_{i}}$
$\displaystyle=\frac{e^{o_{c_{i}}}}{\sum^{|C_{\textit{seen}}|}_{s=1}e^{o_{c_{s}}}}$
(1)

where $|C_{\textit{seen}}|$ is the number of classes that the model has seen,
$l_{c_{i}}\in\\{0,1\\}$ is the one-hot label of class ${c_{i}}$, and
$o(x;\theta,\phi)=[o_{c_{1}},o_{c_{2}},...,o_{c_{|C_{\textit{seen}}|}}]$ is
the set of logit values for input $x$. Given a training sample $x$ of class
$c_{i}$, the gradients on logits ($c_{j}\neq c_{i}$) are given by

$\frac{\partial\mathcal{L}_{ce}(o(x;\theta,\phi))}{\partial
o_{c_{i}}}=p^{c_{i}}-1,\quad\frac{\partial\mathcal{L}_{ce}(o(x;\theta,\phi))}{\partial
o_{c_{j}}}=p^{c_{j}}.$ (2)

From Eq. 2, we see that $x$ gives its true logit $o_{c_{i}}$ a negative
gradient and the other logits $o_{c_{j}}$ positive gradients. As the gradient
update rule for a parameter $w$ is $w=w-lr*\nabla w$, where $lr$ is the
learning rate, the negative gradient ($p^{c_{i}}-1$) results in an increase in
$o_{c_{i}}$ for the true class $c_{i}$ and the positive gradient ($p^{c_{j}}$)
results in a decrease in $o_{c_{j}}$ for each wrong class $c_{j}$. Thus, the
negative gradient encourages the model to output a larger probability for the
true class, and positive gradients help output lower probabilities for the
wrong classes.

However, in CL, as the model has no access to the training data of previous
tasks when it learns a new task incrementally, all gradients on previous
classes are positive during the new task training (there are no negative
gradients) (imbalance of positive and negative gradients). Then the model
tends to output smaller probabilities on previous classes, biasing the
classification towards the new classes.

##  4 Gradient Imbalance (GI) in Replaying

Figure 1: PN rate (rate in the figures) of the CIFAR100 or TinyImageNet
experiments with Adam optimizers. The buffer size is 1000. We choose four
different tasks in each experiment and plot their PN rates as subsequent tasks
are learned. In (C), we report their accuracy. The plots based on the SGD
optimizer are given in Appendix 1, which show the same trend.

A replay-based method actually makes the (positive and negative) gradients
more balanced. There are adjustments of both ups and downs for logits for
previous classes. However, this is insufficient. We discuss two reasons: data
imbalance and CL imbalance. We propose some metrics first.

Metrics: When the model is learning the new task $t$, for a seen class $c_{i}$
in any seen task (including the current task $t$), we define the mean of
positive gradients ($p_{k}^{c_{i}}$) that the logits of all other classes than
$c_{i}$ have and the mean of negative gradients ($p_{k}^{c_{i}}-1$) that class
$c_{i}$’s logit has over all $n^{t,r}$ training samples (including both the
data samples seen so far from the current task $t$ and the corresponding
replay data used) for this task as

$\displaystyle P(t,c_{i})$
$\displaystyle=\frac{\sum_{k=1}^{n^{t,r}}p_{k}^{c_{i}}\cdot\mathbb{I}(y_{k}\neq
c_{i})}{n^{t}}$ (3) $\displaystyle N(t,c_{i})$
$\displaystyle=\frac{\sum_{k=1}^{n^{t,r}}(p_{k}^{c_{i}}-1)\cdot\mathbb{I}(y_{k}=c_{i})}{n^{t}}$

where $y_{k}$ is the label of the $k$th sample of task $t$ and $\mathbb{I}$ is
the indicator function.

At the task level, for classes in any task $t^{\prime}$ ($t^{\prime}\leq t$),
we define the mean of positive and negative gradients that the logits of the
class set $L^{(t^{\prime})}$ receives in training task $t$ respectively as

$\displaystyle P(t,L^{(t^{\prime})})$
$\displaystyle=\frac{\sum_{j=1}^{|C_{\textit{seen}}|}P(t,c_{j})\cdot\mathbb{I}(c_{j}\in
L^{(t^{\prime})})}{\sum_{j=1}^{|C_{\textit{seen}}|}1\cdot\mathbb{I}(c_{j}\in
L^{(t^{\prime})})}$ (4) $\displaystyle N(t,L^{(t^{\prime})})$
$\displaystyle=\frac{\sum_{j=1}^{|C_{\textit{seen}}|}N(t,c_{j})\cdot\mathbb{I}(c_{j}\in
L^{(t^{\prime})})}{\sum_{j=1}^{|C_{\textit{seen}}|}1\cdot\mathbb{I}(c_{j}\in
L^{(t^{\prime})})}$

###  4.1 Gradient Imbalance due to Data Imbalance

In a replay-based method, the memory buffer size is usually fixed and small
and the batch size $N^{\textit{buf}}$ for the buffer batch $X^{\textit{buf}}$
and the batch size $N^{\textit{new}}$ for the new data batch
$X^{\textit{new}}$ are also fixed. As the number of previous classes grows
with more tasks learned, the number of sampled data for each previous class
$\frac{N^{\textit{buf}}}{\bigcup_{r=1}^{t-1}|L^{(r)}|}$ in $X^{\textit{buf}}$
gets smaller, but the number of samples for each new class
$\frac{N^{\textit{new}}}{|L^{(t)}|}$ in $X^{\textit{new}}$ of the new task
remains unchanged. Then we have
$\frac{N^{\textit{new}}}{|L^{(t)}|}>\frac{N^{\textit{buf}}}{\bigcup_{r=1}^{t-1}|L^{(r)}|}$
and the samples from previous classes and the new data classes can become
highly imbalanced. The positive gradients for the previous classes can surpass
their negative gradients (in absolute values). We can empirically verify this
by calculating the positive-negative (PN) gradient rate with $t^{\prime}\leq
t$:

$\textit{PN}(t,t^{\prime})=\frac{P(t,L^{(t^{\prime})})}{N(t,L^{(t^{\prime})})}$
(5)

Empirical Verification: We introduce the setup first.

Datasets and Tasks: We conduct our experiments using two dataset settings: (1)
Split CIFAR100, where we divide the CIFAR100 dataset into 10 tasks with 10
unique classes per task, and (2) Split TinyImagenet, where we divide the
TinyImageNet dataset into 20 tasks with 10 different classes per task. For
online CL, we run each task in one epoch.

Model, Optimizer and Batch Size: We use the full size ResNet-18 to perform our
experiments. To give a more general analysis, we run with both SGD and Adam
optimizers. For learning rate, we follow [3] and set it as 0.1 for the SGD
optimizer. To ensure a good performance, we search and set the learning rate
as 0.001 for the Adam optimizer. For batch size, we also follow [3] and set it
to 10 for $X^{\textit{new}}$ and 10 for the buffer batch $X^{\textit{buf}}$,
which is randomly sampled from the memory buffer. Since ER (Experience Replay)
is a basic replay method, we conduct our experiments using ER [11].

We show the PN results of the classes of a few tasks of the two datasets in
Fig. 1 as more subsequent tasks are learned. We can make the following
observations from Fig. 1. See plots based on SGD in Appendix 1 (see footnote
2).

(1). Sub-figures (A) and (B) in Fig. 1 show that when the model is trained
with the new task $t$, the PN rates for the classes of old tasks get smaller
than $-1$ and tend to decrease. This indicates that positive gradients surpass
negative gradients (in absolute values) for the old classes. The reason is
that the new data classes dominate the training, and the new data give high
positive gradients to classes of the old tasks to ensure the new task data are
not classified to those old classes. As more tasks are learned, the data get
even more imbalanced, resulting in lower negative gradients on the old classes
from the replay data of the old task $t^{\prime}$ and thus the downward trend
of the curves in (A) and (B).

(2). Sub-figure (C) in Fig. 1 shows that the tendency of the test accuracy
performance of the previous classes is similar to the tendency of their PN
rates.

Dealing with Data Imbalance. We propose a sampling strategy that guarantees
the same number of training samples for each class in each training iteration.
Specifically, we fix $N^{\textit{new}}+N^{\textit{buf}}=20$, where
$N^{\textit{new}}$ and $N^{\textit{buf}}$ are the number of samples in
$X^{\textit{new}}$ and $X^{\textit{buf}}$, respectively. We use the ratio
between $|L^{(t)}|$ and $\sum_{r=1}^{t}|L^{(r)}|$ to decide the number of
samples from $X^{\textit{new}}$ (the rest are not used here), i.e.,
$N^{\textit{new}}=max(int(20\cdot\frac{|L^{(t)}|}{\sum_{r=1}^{t}|L^{(r)}|}),1)$
where $int(\cdot)$ returns the nearest integer of a given number and
$N^{\textit{buf}}=20-N^{\textit{new}}$. We use $N^{\textit{buf}}$ to sample
previous tasks’ data in the memory buffer. In this way, the number of training
data for every class is approximately equal.

Figure 2: A-PN rate (rate) for each class ($x$-axis) after the last task is
trained. The classes are ordered by their sequence of appearances in the
original data. The buffer size is 1000 and the optimizer is Adam. (C) shows
the test performances (see the formula in the legion) of all learned classes.
The results based on the SGD optimizer are given in Appendix 3.

###  4.2 Gradient Imbalance due to CL Imbalance

GI still occurs even after the number of samples used in each class is
balanced using the technique above. Let us call the balanced data for training
$X^{\textit{mix}}$, which includes both the replay data from previous tasks
and new data from the current task. We now introduce the second cause of GI,
CL imbalance, which is due to incremental training in CL. To explore CL
imbalance, we define the rate of accumulated positive and negative gradients
(A-PN rate) for a class $c_{i}\in L^{(t^{\prime})}$ from task $t^{\prime}$ to
the current task $t$ ($t^{\prime}\leq t$),

$\small{\textit{A-PN}(t^{\prime},t,c_{i})=\frac{\sum_{r=t^{\prime}}^{t}P(r,c_{i})}{\sum_{r=t^{\prime}}^{t}N(r,c_{i})}}$
(6)

$\textit{A-PN}(t^{\prime},t,c_{i})$ gives the accumulated gradient rate of
each class over the learning process to the current task $t$.

Empirical Verification: We plot the A-PN rate of each class in Fig. 2(A) for
the CIFAR100 dataset and in Fig. 2(B) for the TinyImageNet dataset after the
last task is learned. We make the following observations from Fig. 2 (A), (B),
(C) and some more details in the experiments:

(1). The A-NP rate accumulated is close to -1 or balanced when we jointly
train all classes in the CIFAR100 dataset as a single task (the blue curve and
its moving average yellow curve in Fig. 2 (A)). We also observe balanced
gradients from detailed results (not shown here) for all classes when we learn
the first task (no replay data).

(2). We ignore the accumulative part in Eq. 6 for the time being and focus on
the gradient imbalance in the classes of each task alone. Based on the
detailed results (not given here), we observe that as more tasks are learned,
the GI for the classes of the current task being learned gets worse and worse.
This explains why the classes of the last task in sub-figures (A) and (B) have
the highest GI. Note, the last task has no accumulation as $t^{\prime}=t$. Let
us explain why.

As an example, we consider class $c_{i}$ from the current task $t$. As more
and more tasks are learned, the replay data from previous tasks in
$X^{\textit{mix}}$ will contribute less and less positive gradients to $c_{i}$
because the replay data may have been trained many times in the past and are
well overfitted to their own classes in the past and their probability of
being classified to the new class $c_{i}$ is very small and hence they
contribute very low positive gradients to $c_{i}$. In our experiments, we
observed this and also decreased negative gradient on $c_{i}$. We plot the
average cross-entropy loss of the new data batch and the replay data batch
when the model learns the second task of CIFAR100 in Appendix 2. The figure
shows that the average cross-entropy loss of the new data batch is much higher
than that of the replay data batch in the whole training process. This is
because the system has learned good features from previous tasks that makes
the learning of the new class $c_{i}$ easier resulting in less negative
gradients. But positive gradients drop more significantly.

We now consider the accumulation part (the summations in Eq. 6). The next
section will show that A-PN is used to dynamically adjust the loss for each
class in learning to balance the gradients. By right, when we learn the
current task, we only need to consider the current GI situation to adjust the
loss. However, considering only one case (the current situation) is risky due
to random fluctuations. That is why we consider the impact on a class $c_{i}$
from all subsequent tasks, which gives us a more robust gradient rate for
adjusting the loss in dealing with GI in learning.

(3). Fig. 2(C) shows that higher A-PN rates (imbalanced gradients) result in
lower test accuracy. The fact that the accumulated negative gradient being
greater than the accumulated positive gradients in absolute values (i.e.,
$\textit{A-PN}>-1$) makes the model biased towards the new classes.

##  5 The Proposed GSA Method

Our method GSA consists of two parts: (1) a new optimization framework that
separately optimizes cross-task classification and within-task classification.
(2) A gradient-based self-adaptive loss to alleviate the gradient imbalance
(CL imbalance) in our framework.

###  5.1 Optimizing Within-class Classification and Cross-task Classification

Assume that the system has seen $n-1$ previous classes $(c_{1},...,c_{n-1})$
and there $m-n+1$ new classes ($c_{n},...,c_{m}$) in the current batch
$X^{\textit{new}}$, for a class, say $c_{n}$, in the batch of new data, we
decompose its learning into two parts: (1) learning the decision boundaries
between $c_{n}$ and the other classes in the current batch $X^{\textit{new}}$
and (2) learning the decision boundaries between $c_{n}$ and all previous
$n-1$ classes of old tasks in $X^{\textit{buf}}$. For example, for a sample
$x_{c_{n}}$ of new class $c_{n}$ from $X^{\textit{new}}$, its original cross-
entropy loss is:

$\mathcal{L}_{ce}(x_{c_{n}})=-\log(\frac{e^{o_{c_{n}}}}{\sum^{m}_{s=1}e^{o_{c_{s}}}})$
(7)

where $o_{c_{s}}$ is the logit value of $x_{c_{n}}$ for class $c_{s}$. We
decompose the loss (denoted by $\mathcal{L}_{decom}(x_{c_{n}})$) into two
parts:

$\mathcal{L}_{decom}(x_{c_{n}})=-\log(\frac{e^{o_{c_{n}}}}{\sum^{m}_{s=n}e^{o_{c_{s}}}})-\log(\frac{e^{o_{c_{n}}}}{\sum^{n-1}_{s=1}e^{o_{c_{s}}}+e^{o_{c_{n}}}})$
(8)

The first part of $\mathcal{L}_{decom}$ distinguishes $c_{n}$ and other new
classes with respect to $x_{c_{n}}$. Their gradient rates are similar. The
second part distinguishes $c_{n}$ and the previous classes with respect to
$x_{c_{n}}$. The negative interference from new classes to previous classes is
limited in the second part. For the relationship of loss $\mathcal{L}_{decom}$
and loss $\mathcal{L}_{ce}$, we have the following proposition:

Proposition 1. For a sample $x_{c_{n}}$ from the new class $c_{n}$ in
$X^{\textit{new}}$, the following holds

$\mathcal{L}_{decom}(x_{c_{n}})\geq\mathcal{L}_{ce}(x_{c_{n}})$ (9)

The proof is in Appendix 4. Similarly, for a sampled data point $x_{c_{j}}$ of
class $c_{j}(j<n)$ from the buffer data $X^{\textit{buf}}$, we define its new
loss as:

$\mathcal{L}_{decom}(x_{c_{j}})=-\log(\frac{e^{o_{c_{j}}}}{\sum^{n-1}_{s=1}e^{o_{c_{s}}}})-\log(\frac{e^{o_{c_{j}}}}{\sum^{m}_{s=n}e^{o_{c_{s}}}+e^{o_{c_{j}}}})$
(10)

In this loss, we put the previous classes (similar gradient rates) in the
first term to maintain the established decision boundaries between previous
classes and put the logit of the true class and the logits of the classes in
the current task in the second term to establish the decision boundaries
between class $c_{j}$ and the new classes. The negative transfer from the new
classes to previous classes is limited to the second term. We prove that
$\mathcal{L}_{decom}(x_{c_{j}})$ is the upper bound of
$\mathcal{L}_{ce}(x_{c_{j}})$ in Appendix 4.

Considering data balance. To estimate the first term in Eq. 8 for class
$c_{n}$, we use $X^{\textit{new}}$ as the data of new classes as
$X^{\textit{new}}$ is usually class-balanced.111If this is not the case, we
can sample some data of the same class from the memory buffer or directly
duplicate the samples to guarantee it. For a similar reason, we use
$X^{\textit{buf}}$ to estimate the first term in Eq. 10 for the previous class
$c_{j}$. However, for the second term in Eq. 8, using $X^{\textit{new}}$
causes prediction bias towards new classes as $X^{\textit{new}}$ lacks samples
of previous classes. Also, using $X^{\textit{buf}}$ to estimate the second
term in Eq. 10 is not appropriate as the number of new classes samples in the
buffer increases gradually and for a long time, the new classes are the minor
classes in the buffer. To solve both problems, we use the sampling strategy
proposed in the last section under “Dealing with Data Imbalance.” Then we
obtain a mixed set of new data and old data sampled from the memory buffer,
called $X^{\textit{mix}}$. $X^{\textit{mix}}$ can be regarded as samples from
the uniform joint distribution for all classes, including both the new and old
classes. Then our final optimization goal is:

$\displaystyle\min\limits_{F\in\mathcal{F}}-\mathop{\mathbb{E}}\limits_{x\sim
X^{\textit{new}}}\log(\frac{e^{o_{y_{x}}}}{\sum^{m}_{s=n}e^{o_{c_{s}}}})-\mathop{\mathbb{E}}\limits_{x\sim
X^{\textit{buf}}}\log(\frac{e^{o_{y_{x}}}}{\sum^{n-1}_{s=1}e^{o_{c_{s}}}})$
(11) $\displaystyle-\mathop{\mathbb{E}}\limits_{x\sim
X^{\textit{mix}}}({\mathbb{I}(y_{x}\cap
S(X^{\textit{new}})\neq\emptyset)}\cdot\log(\frac{e^{o_{y_{x}}}}{\sum^{n-1}_{s=1}e^{o_{c_{s}}}+e^{o_{y_{x}}}})+$
$\displaystyle{\mathbb{I}(y_{x}\cap
S(X^{\textit{new}})=\emptyset)}\cdot\log(\frac{e^{o_{y_{x}}}}{\sum^{m}_{s=n}e^{o_{c_{s}}}+e^{o_{y_{x}}}}))$

where $y_{x}$ is the label of sample $x$ and $S(X^{\textit{new}})$ is the set
of classes appeared in the $X^{\textit{new}}$ and $\mathbb{I}$ is the
indicator function and $\mathcal{F}$ is the function space. The first term and
second terms focus on the within-task classification goal and the last term
focus on the cross-task classification goal.

###  5.2 Self-Adaptive Loss for CL Imbalance

Based on the analysis in Section 4, we know that even though the data balance
is guaranteed, GI still exists as the gradient rate is different for different
classes, which is caused by CL imbalance. We note that CL imbalance (measured
by the gradient rate of each class) not only occurs between new and previous
classes but also occurs among previous classes. To mitigate it, we propose
GSA-CE (Gradient-based Self-Adaptive CE) loss $\mathcal{L}_{\textit{GSA-CE}}$.
The loss function is dynamically adjusted based on the gradient rates, which
change as more tasks and batches are learned. Let
$\\{x_{k},y_{k}\\}_{k=1}^{N^{\textit{mix}}}$ be the samples in
$X^{\textit{mix}}$, their $\mathcal{L}_{\textit{GSA-CE}}$ is:

$\displaystyle\mathcal{L}_{\textit{GSA-
CE}}(\\{x_{k},y_{k}\\}_{k=1}^{N^{\textit{mix}}}=-\frac{1}{N^{\textit{mix}}}\sum_{k=1}^{N^{\textit{mix}}}w_{y_{k}}\cdot({\mathbb{I}(y_{k}\cap
S(X^{\textit{new}})\neq\emptyset)}\cdot$ (12)
$\displaystyle\log(\frac{e^{o_{y_{k}}}}{\sum^{n-1}_{s=1}e^{o_{c_{s}}}+e^{o_{y_{k}}}})+{\mathbb{I}(y_{k}\cap
S(X^{\textit{new}})=\emptyset)}\cdot$
$\displaystyle\log(\frac{e^{o_{y_{k}}}}{\sum^{m}_{s=n}v_{c_{s}}e^{o_{c_{s}}}+e^{o_{y_{k}}}}))$

where

$\small w_{y_{k}}=\frac{2}{1+e^{\textit{A-PN}(t^{\prime},t,y_{k})}},\quad
v_{c_{s}}=\frac{1}{-\textit{PN}(t,c_{s})}$ (13)

where $t^{\prime}$ is the task in which class $c_{s}$ appears. When the model
is training task $t$, we incrementally update the
$\textit{A-PN}(t^{\prime},t,y_{k})$ and $\textit{PN}(t,y_{k})$ by adding the
new gradients into the current sums of all previous gradients.
$\mathcal{L}_{\textit{GSA-CE}}$ has three advantages:

(1). When the new task arrives, we do not need to assume that we have the
entire training set for calculating some statistics or an exemplar set like
CCIL [33]. Our method is thus suitable for online CL.

(2). The model sets the values of $w_{y_{k}}$ and $v_{c_{s}}$ based on A-PN
rates and PN rates automatically as they vary dynamically with the current
gradient imbalance situation.

(3). A-PN and PN rates are class-based metrics that reflect the discriminative
power of the model for all seen classes. Our loss then gives class-based
weights for samples from different classes. The new loss thus automatically
adjust the loss weight $w_{y_{k}}$ to balance the accumulated negative
gradient and the accumulated positive gradient for each class. We justify the
new loss in Appendix 5.

(4) No new hyperparameter is introduced in the process.

Similarly, as the gradient rate imbalance occurs also within previous classes,
for samples $\\{x_{k},y_{k}\\}_{k=1}^{N^{\textit{buf}}}$ from
$X^{\textit{buf}}$, their $\mathcal{L}_{\textit{GSA-CE}}$ loss is:

$\displaystyle\mathcal{L}_{\textit{GSA-
CE}}(\\{x_{k},y_{k}\\}_{k=1}^{N^{\textit{buf}}})=-\frac{\sum_{k=1}^{N^{\textit{buf}}}w_{y_{k}}\log(\frac{e^{o_{y_{k}}}}{\sum^{n-1}_{s=1}e^{o_{c_{s}}}})}{N^{\textit{buf}}}$
(14)

Considering the gradient imbalance, we replace the second term and the third
term in Eq. 11 with our $\mathcal{L}_{\textit{GSA-CE}}$ losses Eq. 14 and Eq.
12 in the expectation form respectively. The final optimization objective in
Eq. 11 is turned into:

$\displaystyle\min\limits_{F\in\mathcal{F}}-\mathop{\mathbb{E}}\limits_{x\sim
X^{\textit{new}}}\log(\frac{e^{o_{y_{x}}}}{\sum^{m}_{s=n}e^{o_{c_{s}}}})-\mathop{\mathbb{E}}\limits_{x\sim
X^{\textit{buf}}}w_{y_{x}}\cdot\log(\frac{e^{o_{y_{x}}}}{\sum^{n-1}_{s=1}e^{o_{c_{s}}}})$
(15) $\displaystyle-\mathop{\mathbb{E}}\limits_{x\sim
X^{\textit{mix}}}w_{y_{x}}\cdot({\mathbb{I}(y_{x}\cap
S(X^{\textit{new}})\neq\emptyset)}\cdot\log(\frac{e^{o_{y_{x}}}}{\sum^{n-1}_{s=1}e^{o_{c_{s}}}+e^{o_{y_{x}}}})$
$\displaystyle+{\mathbb{I}(y_{x}\cap
S(X^{\textit{new}})=\emptyset)}\cdot\log(\frac{e^{o_{y_{x}}}}{\sum^{m}_{s=n}v_{c_{s}}e^{o_{c_{s}}}+e^{o_{y_{x}}}}))$

##  6 Experiments

Method |  | MNIST (5 tasks) |  |  | CIFAR10 (5 tasks) |  |  | CIFAR100 (10 tasks) |  |  | TinyImageNet (100 tasks) |
---|---|---|---|---|---|---|---|---|---|---|---|---
M |  M=0.1k |  M=0.5k |  M=1k |  M=0.2k |  M=0.5k |  M=1k |  M=1k |  M=2k |  M=5k |  M=2k |  M=4k |  M=10k
AGEM [10] |  56.9$\pm$5.2 |  57.7$\pm$8.8 |  61.6$\pm$3.2 |  22.7$\pm$1.8 |  22.7$\pm$1.9 |  22.6$\pm$0.7 |  5.8$\pm$0.2 |  6.5$\pm$0.3 |  5.8$\pm$0.2 |  0.9$\pm$0.1 |  2.1$\pm$0.1 |  3.9$\pm$0.2
ER [11] |  78.7$\pm$0.4 |  88.0$\pm$0.2 |  90.3$\pm$0.1 |  49.7$\pm$0.6 |  55.2$\pm$0.6 |  59.3$\pm$0.2 |  15.7$\pm$0.3 |  22.4$\pm$0.5 |  29.5$\pm$0.9 |  4.7$\pm$0.5 |  10.1$\pm$0.7 |  11.7$\pm$0.2
MIR [3] |  79.0$\pm$0.5 |  88.3$\pm$0.1 |  91.3$\pm$1.9 |  37.3$\pm$0.3 |  40.0$\pm$0.6 |  41.0$\pm$0.6 |  15.7$\pm$0.2 |  19.1$\pm$0.1 |  24.1$\pm$0.2 |  6.1$\pm$0.5 |  11.7$\pm$0.2 |  13.5$\pm$0.2
GSS [4] |  70.4$\pm$1.5 |  80.7$\pm$5.8 |  87.5$\pm$5.9 |  26.9$\pm$1.2 |  30.7$\pm$1.3 |  40.1$\pm$1.4 |  11.1$\pm$0.2 |  13.3$\pm$0.5 |  17.4$\pm$0.1 |  3.3$\pm$0.5 |  10.0$\pm$0.2 |  10.5$\pm$0.2
ASER [40] |  61.6$\pm$2.1 |  71.0$\pm$0.6 |  82.1$\pm$5.9 |  27.8$\pm$1.0 |  36.2$\pm$1.2 |  44.7$\pm$1.2 |  16.4$\pm$0.3 |  12.2$\pm$1.9 |  27.1$\pm$0.3 |  5.3$\pm$0.3 |  8.2$\pm$0.2 |  10.3$\pm$0.4
Rainbow-rt [5] |  89.1$\pm$0.3 |  92.1$\pm$0.1 |  95.0$\pm$0.3 |  45.2$\pm$0.4 |  50.6$\pm$0.3 |  51.5$\pm$0.8 |  15.4$\pm$0.5 |  15.9$\pm$0.3 |  20.4$\pm$0.3 |  6.6$\pm$0.3 |  10.1$\pm$0.3 |  13.1$\pm$0.5
ER-AML [7] |  76.5$\pm$0.1 |  86.0$\pm$0.2 |  91.5$\pm$0.1 |  40.5$\pm$0.7 |  48.7$\pm$0.5 |  50.1$\pm$0.4 |  16.1$\pm$0.4 |  17.6$\pm$0.5 |  22.6$\pm$0.1 |  5.4$\pm$0.2 |  7.1$\pm$0.5 |  10.1$\pm$0.4
GDumb [35] |  81.2$\pm$0.5 |  91.0$\pm$0.2 |  94.5$\pm$0.1 |  35.9$\pm$1.1 |  50.7$\pm$0.7 |  63.5$\pm$0.5 |  17.1$\pm$0.4 |  25.1$\pm$0.2 |  38.6$\pm$0.5 |  12.6$\pm$0.1 |  12.7$\pm$0.3 |  15.7$\pm$0.2
SCR [31] |  86.2$\pm$0.5 |  92.8$\pm$0.3 |  94.6$\pm$0.1 |  47.2$\pm$1.7 |  58.2$\pm$0.5 |  64.1$\pm$1.2 |  26.5$\pm$0.2 |  31.6$\pm$0.5 |  36.5$\pm$0.2 |  10.6$\pm$1.1 |  17.2$\pm$0.1 |  20.4$\pm$1.1
DER++ [6] |  74.4$\pm$1.1 |  91.5$\pm$0.2 |  92.1$\pm$0.2 |  44.2$\pm$1.1 |  47.9$\pm$1.5 |  54.7$\pm$2.2 |  15.3$\pm$0.2 |  19.7$\pm$1.5 |  27.0$\pm$0.7 |  4.5$\pm$0.3 |  10.1$\pm$0.3 |  17.6$\pm$0.5
IL2A [55] |  90.2$\pm$0.1 |  92.7$\pm$0.1 |  93.9$\pm$0.1 |  54.7$\pm$0.5 |  56.0$\pm$0.4 |  58.2$\pm$1.2 |  18.2$\pm$1.2 |  19.7$\pm$0.5 |  22.4$\pm$0.2 |  5.5$\pm$0.7 |  8.1$\pm$1.2 |  11.6$\pm$0.4
Co2L [9] |  83.1$\pm$0.1 |  91.5$\pm$0.1 |  94.7$\pm$0.1 |  42.1$\pm$1.2 |  51.0$\pm$0.7 |  58.8$\pm$0.4 |  17.1$\pm$0.4 |  24.2$\pm$0.2 |  32.2$\pm$0.5 |  10.1$\pm$0.2 |  15.8$\pm$0.4 |  22.5$\pm$1.2
LUCIR [19] |  73.2$\pm$0.1 |  87.2$\pm$0.2 |  90.1$\pm$0.1 |  27.9$\pm$1.2 |  33.5$\pm$0.7 |  37.8$\pm$0.5 |  8.6$\pm$1.3 |  19.5$\pm$0.7 |  16.9$\pm$0.5 |  7.6$\pm$0.5 |  9.6$\pm$0.7 |  12.5$\pm$0.7
CCIL [33] |  86.4$\pm$0.1 |  92.8$\pm$0.2 |  94.0$\pm$0.1 |  50.5$\pm$0.2 |  55.3$\pm$0.5 |  59.8$\pm$0.3 |  18.5$\pm$0.3 |  19.1$\pm$0.4 |  20.5$\pm$0.3 |  5.6$\pm$0.9 |  7.0$\pm$0.5 |  15.2$\pm$0.5
BiC [46] |  90.4$\pm$0.1 |  93.0$\pm$0.2 |  94.8$\pm$0.1 |  48.2$\pm$0.7 |  57.5$\pm$1.4 |  63.8$\pm$0.2 |  21.2$\pm$0.3 |  36.1$\pm$1.3 |  42.5$\pm$1.2 |  10.2$\pm$0.9 |  18.9$\pm$0.3 |  25.2$\pm$0.6
SSIL [2] |  88.2$\pm$0.1 |  93.0$\pm$0.2 |  95.1$\pm$0.1 |  49.5$\pm$0.2 |  59.2$\pm$0.4 |  64.0$\pm$0.5 |  26.0$\pm$0.1 |  33.1$\pm$0.5 |  39.5$\pm$0.4 |  9.6$\pm$0.7 |  15.2$\pm$1.5 |  21.1$\pm$0.1
GSA |  91.4$\pm$0.1 |  93.2$\pm$0.1 |  96.5$\pm$0.1 |  58.0$\pm$0.4 |  64.6$\pm$0.2 |  69.1$\pm$0.3 |  31.4$\pm$0.2 |  39.7$\pm$0.6 |  49.7$\pm$0.2 |  18.4$\pm$0.4 |  26.0$\pm$0.2 |  33.2$\pm$0.4

Table 1: Accuracy on the four experiment datasets with different memory buffer
sizes M. All values are averages of 15 runs.

Method |  | MNIST (5 tasks) |  |  | CIFAR10 (5 tasks) |  |  | CIFAR100 (10 tasks) |  |  | TinyImageNet (100 tasks) |
---|---|---|---|---|---|---|---|---|---|---|---|---
M |  M=0.1k |  M=0.5k |  M=1k |  M=0.2k |  M=0.5k |  M=1k |  M=1k |  M=2k |  M=5k |  M=2k |  M=4k |  M=10k
AGEM |  32.5$\pm$5.9 |  30.1$\pm$4.2 |  32.0$\pm$2.9 |  36.1$\pm$3.8 |  43.2$\pm$4.3 |  48.1$\pm$3.4 |  78.6$\pm$2.1 |  77.5$\pm$1.3 |  78.3$\pm$1.2 |  73.9$\pm$0.2 |  78.9$\pm$0.2 |  74.1$\pm$0.3
ER |  22.7$\pm$0.5 |  9.7$\pm$0.4 |  6.7$\pm$0.5 |  42.0$\pm$0.3 |  26.7$\pm$0.7 |  20.7$\pm$0.7 |  65.1$\pm$1.3 |  59.3$\pm$0.9 |  60.0$\pm$1.6 |  68.2$\pm$2.8 |  66.2$\pm$0.8 |  67.2$\pm$0.2
MIR [3] |  22.3$\pm$0.5 |  9.0$\pm$0.5 |  5.7$\pm$0.9 |  40.0$\pm$1.6 |  25.9$\pm$0.7 |  24.5$\pm$0.5 |  24.5$\pm$0.3 |  21.4$\pm$0.3 |  21.0$\pm$0.1 |  61.1$\pm$3.2 |  60.9$\pm$0.3 |  59.5$\pm$0.3
GSS |  26.1$\pm$2.2 |  17.8$\pm$5.22 |  10.5$\pm$6.7 |  75.5$\pm$1.5 |  65.9$\pm$1.6 |  54.9$\pm$2.0 |  73.4$\pm$4.2 |  69.3$\pm$3.1 |  70.9$\pm$2.9 |  72.8$\pm$1.2 |  72.6$\pm$0.4 |  71.5$\pm$0.2
ASER [40] |  33.8$\pm$1.1 |  24.8$\pm$0.5 |  13.8$\pm$0.4 |  71.1$\pm$1.8 |  59.1$\pm$1.5 |  50.4$\pm$1.5 |  25.0$\pm$0.2 |  12.2$\pm$1.9 |  13.2$\pm$0.1 |  65.7$\pm$0.7 |  64.2$\pm$0.2 |  62.2$\pm$0.1
Rainbow-rt |  10.1$\pm$0.1 |  4.7$\pm$0.4 |  2.4$\pm$0.5 |  20.4$\pm$0.3 |  18.1$\pm$0.4 |  15.3$\pm$0.8 |  25.5$\pm$0.5 |  19.3$\pm$0.4 |  13.3$\pm$0.5 |  25.5$\pm$0.6 |  23.2$\pm$0.3 |  20.0$\pm$0.4
ER-AML |  23.1$\pm$0.1 |  9.7$\pm$0.4 |  6.4$\pm$0.5 |  50.9$\pm$0.3 |  40.1$\pm$0.5 |  34.2$\pm$0.8 |  51.5$\pm$0.8 |  49.2$\pm$0.5 |  38.7$\pm$0.6 |  47.4$\pm$0.5 |  43.2$\pm$0.3 |  41.0$\pm$0.5
GDumb |  10.3$\pm$0.1 |  6.2$\pm$0.1 |  4.8$\pm$0.2 |  26.5$\pm$0.5 |  24.5$\pm$0.2 |  18.9$\pm$0.4 |  16.7$\pm$0.5 |  17.6$\pm$0.2 |  16.8$\pm$0.4 |  15.9$\pm$0.5 |  14.6$\pm$0.3 |  11.7$\pm$0.2
SCR |  10.7$\pm$0.1 |  4.7$\pm$0.1 |  4.0$\pm$0.2 |  41.3$\pm$0.1 |  31.5$\pm$0.2 |  24.7$\pm$0.4 |  17.5$\pm$0.2 |  11.6$\pm$0.5 |  5.6$\pm$0.4 |  19.4$\pm$0.3 |  15.4$\pm$0.3 |  14.9$\pm$0.7
DER++ [6] |  25.0$\pm$0.3 |  7.3$\pm$0.3 |  6.6$\pm$1.2 |  30.1$\pm$0.8 |  31.8$\pm$2.5 |  18.7$\pm$3.4 |  43.4$\pm$0.2 |  44.0$\pm$1.9 |  25.8$\pm$3.5 |  67.2$\pm$1.7 |  63.6$\pm$0.3 |  55.2$\pm$0.7
IL2A [55] |  8.7$\pm$0.1 |  7.2$\pm$0.1 |  4.1$\pm$0.1 |  36.0$\pm$0.2 |  32.1$\pm$0.4 |  29.1$\pm$0.4 |  24.6$\pm$0.6 |  12.5$\pm$0.7 |  20.0$\pm$0.5 |  65.5$\pm$0.7 |  60.1$\pm$0.5 |  57.6$\pm$1.1
Co2L [9] |  14.7$\pm$0.2 |  7.1$\pm$0.1 |  3.1$\pm$0.1 |  32.0$\pm$0.1 |  21.0$\pm$0.3 |  16.9$\pm$0.2 |  16.9$\pm$0.4 |  16.6$\pm$0.6 |  9.9$\pm$0.7 |  60.5$\pm$0.5 |  52.5$\pm$0.9 |  42.5$\pm$0.8
LUCIR |  20.3$\pm$0.1 |  8.5$\pm$0.1 |  7.8$\pm$0.1 |  63.5$\pm$0.5 |  55.3$\pm$0.2 |  46.5$\pm$1.2 |  60.0$\pm$0.1 |  47.5$\pm$0.9 |  44.3$\pm$0.7 |  46.4$\pm$0.7 |  42.2$\pm$0.9 |  37.6$\pm$0.7
CCIL |  14.1$\pm$0.1 |  7.7$\pm$0.1 |  4.8$\pm$0.1 |  18.6$\pm$0.1 |  16.5$\pm$0.4 |  12.5$\pm$0.8 |  16.7$\pm$0.5 |  16.1$\pm$0.3 |  17.5$\pm$0.2 |  59.4$\pm$0.3 |  56.2$\pm$1.3 |  48.9$\pm$0.6
BiC |  11.1$\pm$0.1 |  3.7$\pm$0.1 |  2.5$\pm$0.1 |  35.4$\pm$0.7 |  25.3$\pm$0.4 |  14.5$\pm$0.7 |  40.2$\pm$0.4 |  30.9$\pm$0.7 |  18.7$\pm$0.5 |  43.5$\pm$0.5 |  32.9$\pm$0.5 |  24.9$\pm$0.4
SSIL |  11.3$\pm$0.1 |  2.7$\pm$0.1 |  2.8$\pm$0.1 |  36.0$\pm$0.7 |  29.6$\pm$0.4 |  13.5$\pm$0.4 |  40.1$\pm$0.5 |  33.9$\pm$1.2 |  21.7$\pm$0.8 |  44.4$\pm$0.7 |  36.6$\pm$0.7 |  29.0$\pm$0.7
GSA |  8.1$\pm$0.1 |  2.5$\pm$0.1 |  1.4$\pm$0.1 |  23.5$\pm$0.2 |  12.6$\pm$0.4 |  10.0$\pm$0.3 |  33.2$\pm$0.6 |  22.8$\pm$0.4 |  8.7$\pm$0.3 |  35.5$\pm$0.3 |  25.8$\pm$0.4 |  16.9$\pm$0.6

Table 2: Average forgetting rate. All numbers are the averages of 15 runs.

Dataset | no new loss | no previous loss | no cross loss |  no $X^{\textit{mix}}$ | no CL imbalance
---|---|---|---|---|---
CIFAR100 |  30.1$\pm$0.6 |  12.7$\pm$0.2 |  29.9$\pm$0.4 |  30.1$\pm$0.1 |  29.7$\pm$0.3
TinyImageNet |  17.5$\pm$0.2 |  4.7$\pm$0.5 |  17.2$\pm$0.4 |  17.5$\pm$0.5 |  16.9$\pm$0.2

Table 3: Ablation accuracy - average of 5 runs. Memory size M is 1k for
CIFAR100 and 2k for TinyImageNet.

Evaluation data. Four image classification datasets are used. 1) MNIST [27]
has 10 classes with 60,000/10,000 training/test examples. We created 5
disjoint tasks with 2 classes per task. 2) CIFAR10 [25] has 10 classes with
50,000/10,000 training/test samples. We created 5 disjoint tasks with 2
classes per task. 3) CIFAR100 [25] has 100 classes with 50,000/10,000
training/test samples. 10 disjoint tasks are created with 10 classes per task.
4) TinyImageNet [26] has 200 classes. We created 100 disjoint tasks with 2
classes per task. Each class has 500 training examples and 50 test examples.

Compared Baselines. GSA222Code and Appendix: https://github.com/gydpku/GSA is
compared with 16 recent baselines: 9 online CL baselines, AGEM, ER, MIR, GSS,
ASER, RainBow-rt, ER-AML, GDumb, and SCR and 7 offline CL baselines, DER++,
IL2A, Co2L, LUCIR, CCIL, BiC and SSIL, as they deal with data/class imbalance
and can be run in the online CL mode, i.e., training in one epoch without
requiring the full data of each task to be available when the task arrives.
The citations of these systems are given in Table 1 associated with their
results.  Note that Rainbow [5] is not a standard online CL method as it
retrains all replay-data with data augmentation for 256 epochs after each
task, which is not suitable for online CL because it cannot be used for any-
time inference. We thus removed this retraining operation and denote it as
Rainbow-rt.

###  6.1 Architecture, Data Augmentation, Training Details and Evaluation
Protocol

Architecture. For MNIST, GSA and baselines employ a fully-connected network
with two hidden layers as the feature extractor $h_{\theta}$, each comprising
of 400 ReLU units. A linear layer of size [400, 10] is used as the classifier
$f_{\phi}$. For CIFAR10, CIFAR100, and TinyImageNet, we follow [6] and use
ResNet18 (not pre-trained) as the feature extractor $h_{\theta}$ with around
11 million trainable parameters for GSA and all baselines. Denoting
$C_{\textit{num}}$ as the number of all classes, we employ a linear layer of
size [$dim_{h}$, $C_{\textit{num}}$] as the classifier $f_{\phi}$. For an
input $x$, we use $F(x)$ to compute $\mathcal{L}_{ce}$.

Data Augmentation. To learn better features, we apply two data augmentations
to each image in $X^{\textit{new}}$ and $X^{\textit{mix}}$: random-resized-
crop and random-gray-scale. For fair comparison, the same data augmentations
are applied to all baselines ($X^{\textit{new}}$ and $X^{\textit{buf}}$) to
improve their performance, which results in an average of about 2% of
improvement for all methods with no drop in performance for any baseline.

Training and hyperparameter settings. Like ER and many other online CL
systems, GSA uses reservoir sampling for memory update. We follow [16] and
train GSA with the Adam optimizer. We set the learning rate as 0.001 and fix
the weight decay as 0.0001 for all settings. We follow [16] and set the batch
size of $X^{\textit{new}}$ as 10 for all methods and $X^{\textit{mix}}$ as 64
for our GSA and $X^{\textit{buf}}$ also as 64 for baselines. We list other
hyper-parameters in Appendix 6.

We use the official codes of baselines. Their source links and default hyper-
parameters are listed in Appendix 6. We run all methods with one epoch for
each task.

Evaluation protocol. We first learn all tasks from the data stream for each
dataset, and then test the final model using the test data of all tasks. We
report the average accuracy of 15 random runs.

###  6.2 Results Analysis and Ablation Experiments

The results in Table 1 show that the best online CL baseline is SCR and the
best class-imbalance baseline is BiC. Although BiC and SSIL were not
originally designed for online CL, they perform well for online CL, better
than all online CL baselines, which do not deal with data imbalance. This
indicates that the data-balanced approaches in BiC and SSIL help. Our method
GSA consistently outperforms them by very large margins as we also deal with
the proposed CL imbalance. Note that OCM is not compared in the tables as it
is not a competitor of GSA. GSA is actually complementary to OCM. For example,
with the largest memory size for each dataset, OCM+GSA gives 96.5% on MNIST,
77.5% on CIFAR10, 53.7% on CIFAR100, and 35.7% on TinyImageNet and also
outperforms OCM [16].

Forgetting rate. Table 2 shows that our GSA has substantially lower forgetting
rates than baselines except for GDumb, Rainbow-rt, and SCR (in two datasets),
but GDumb, Rainbow-rt, and SCR’s accuracy values are substantially lower than
that of GSA (see Table 1). The calculation method of the average forgetting
rate and Training times are given in Appendix 7.

Ablation. We conducted ablation experiments on GSA using two datasets,
CIFAR100 and TinyImageNet. The results are given in Table 3.

(1). Ablation study of training loss in GSA. In the experiments “no new loss”,
”no previous loss”, and ”no cross loss”, we do not consider the first term
(establishing the boundaries between new classes), the second term
(maintaining the decision boundaries between previous classes), and the last
term (establishing the cross-task boundaries) in Eq. 15 respectively. Table 3
shows that their performances are all poorer than GSA (Table 1). This is
because in the first case, the model forgets the knowledge of previous
classes, which causes a drastic drop in performance. The second experiment
shows that our method benefits from considering all the new data in
$X^{\textit{new}}$. The third experiment shows that establishing cross-task
boundaries is an important problem for CIL and can improve the overall
performance further.

(2). Ablation study of data imbalance and CL imbalance. In the experiments “no
balanced sampling,” we replace the sampling strategy for $X^{\textit{mix}}$
with random sampling from the memory buffer $X^{\textit{buf}}$. In the
experiment “no CL imbalance,” we do not use the $\mathcal{L}_{\textit{GSA-
CE}}$ losses proposed in Eq. 12 and Eq. 14 and introduced into the final
optimization objective in Eq. 15, but replace them with the simple cross-
entropy loss. Table 3 shows that all these incomplete GSA systems are poorer
than the full GSA in Table 1. The poorer performance of “no balanced-sampling”
is mainly due to data imbalance between new classes and previous classes in
the buffer. Specifically, in the initial training process of the new task, the
stored training data of the new task classes are usually fewer than that of
the previous classes in the buffer and the cross-task class boundaries are not
established well. The poorer performance of “no CL imbalance” is because the
method does not consider different gradient rates across different classes,
which again results in poorer performances.

GSA-CE loss mitigates the gradient imbalance. Fig. 3 shows that the GSA-CE
loss clearly reduces the A-PN gradient rate imbalance than the method without
it by dropping the A-NP rate that is bigger than -1 and improving the A-NP
rate that is smaller than -1. And it makes the curve of A-PN rate more similar
to the joint training curve.

Figure 3: The gradient rate after the model has learned all classes of
CIFAR100. A-PN is the gradient rate. The memory buffer size is 1000. Joint
training means the model learns all classes together as one task in one epoch.

Long-tail online CL experiments. See Appendix 8, where we will show that our
method also outperforms the baselines.

##  7 Conclusion

This paper discussed the challenge of cross-task class discrimination (CTCD)
and showed that the replay approach partially deals with the problem. However,
the replay approach has a major bias in training, which is manifested by
gradient imbalance on the logits and significantly limits the online CL
performance and the ability to solve the CTCD problem. The paper then analyzed
gradient imbalance from two perspectives: data imbalance introduced by replay
and CL imbalance due to CL itself. After that, it proposed a new learning
strategy and a new self-adaptive loss function to deal with the problems.
Empirical evaluation demonstrated that the new approach GSA improves the
baselines by large margins.

## References

  * [1] Davide Abati, Jakub Tomczak, Tijmen Blankevoort, Simone Calderara, Rita Cucchiara, and Babak Ehteshami Bejnordi.  Conditional channel gated networks for task-aware continual learning.  In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, pages 3931–3940, 2020.
  * [2] Hongjoon Ahn, Jihwan Kwak, Subin Lim, Hyeonsu Bang, Hyojun Kim, and Taesup Moon.  Ss-il: Separated softmax for incremental learning.  In Proceedings of the IEEE/CVF International Conference on Computer Vision, pages 844–853, 2021.
  * [3] Rahaf Aljundi, Lucas Caccia, Eugene Belilovsky, Massimo Caccia, Min Lin, Laurent Charlin, and Tinne Tuytelaars.  Online continual learning with maximally interfered retrieval.  arXiv preprint arXiv:1908.04742, 2019.
  * [4] Rahaf Aljundi, Min Lin, Baptiste Goujaud, and Yoshua Bengio.  Gradient based sample selection for online continual learning.  arXiv preprint arXiv:1903.08671, 2019.
  * [5] Jihwan Bang, Heesu Kim, YoungJoon Yoo, Jung-Woo Ha, and Jonghyun Choi.  Rainbow memory: Continual learning with a memory of diverse samples.  In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, pages 8218–8227, 2021.
  * [6] Pietro Buzzega, Matteo Boschini, Angelo Porrello, Davide Abati, and Simone Calderara.  Dark experience for general continual learning: a strong, simple baseline.  arXiv preprint arXiv:2004.07211, 2020.
  * [7] Lucas Caccia, Rahaf Aljundi, Nader Asadi, Tinne Tuytelaars, Joelle Pineau, and Eugene Belilovsky.  Reducing representation drift in online continual learning.  arXiv preprint arXiv:2104.05025, 2021.
  * [8] Francisco M Castro, Manuel J Marín-Jiménez, Nicolás Guil, Cordelia Schmid, and Karteek Alahari.  End-to-end incremental learning.  In ECCV, 2018.
  * [9] Hyuntak Cha, Jaeho Lee, and Jinwoo Shin.  Co2l: Contrastive continual learning.  In Proceedings of the IEEE/CVF International Conference on Computer Vision (ICCV), pages 9516–9525, October 2021.
  * [10] Arslan Chaudhry, Marc’Aurelio Ranzato, Marcus Rohrbach, and Mohamed Elhoseiny.  Efficient lifelong learning with a-gem.  arXiv preprint arXiv:1812.00420, 2018.
  * [11] Arslan Chaudhry, Marcus Rohrbach, Mohamed Elhoseiny, Thalaiyasingam Ajanthan, Puneet K Dokania, Philip HS Torr, and M Ranzato.  Continual learning with tiny episodic memories.  In ICML-2019, 2019.
  * [12] Arslan Chaudhry, Marcus Rohrbach, Mohamed Elhoseiny, Thalaiyasingam Ajanthan, Puneet K Dokania, Philip HS Torr, and Marc’Aurelio Ranzato.  On tiny episodic memories in continual learning.  arXiv preprint arXiv:1902.10486, 2019.
  * [13] Cyprien de Masson d’Autume, Sebastian Ruder, Lingpeng Kong, and Dani Yogatama.  Episodic memory in lifelong language learning.  In NeurIPS, 2019.
  * [14] Mehrdad Farajtabar, Navid Azizan, Alex Mott, and Ang Li.  Orthogonal gradient descent for continual learning.  In International Conference on Artificial Intelligence and Statistics, pages 3762–3773. PMLR, 2020.
  * [15] Yiduo Guo, Wenpeng Hu, Dongyan Zhao, and Bing Liu.  Adaptive orthogonal projection for batch and online continual learning.  In Proceedings of AAAI-2021, 2022.
  * [16] Yiduo Guo, Bing Liu, and Dongyan Zhao.  Online continual learning through mutual information maximization.  In International Conference on Machine Learning, pages 8109–8126. PMLR, 2022.
  * [17] Yunhui Guo, Mingrui Liu, Tianbao Yang, and Tajana Rosing.  Improved schemes for episodic memory-based lifelong learning.  Advances in Neural Information Processing Systems, 33:1023–1035, 2020.
  * [18] Tyler L Hayes, Kushal Kafle, Robik Shrestha, Manoj Acharya, and Christopher Kanan.  Remind your neural network to prevent catastrophic forgetting.  arXiv preprint arXiv:1910.02509, 2019.
  * [19] Saihui Hou, Xinyu Pan, Chen Change Loy, Zilei Wang, and Dahua Lin.  Learning a unified classifier incrementally via rebalancing.  In CVPR, pages 831–839, 2019.
  * [20] Steven C. Y. Hung, Cheng-Hao Tu, Cheng-En Wu, Chien-Hung Chen, Yi-Ming Chan, and Chu-Song Chen.  Compacting, picking and growing for unforgetting continual learning.  In NeurIPS, 2019.
  * [21] Ronald Kemker and Christopher Kanan.  FearNet: Brain-Inspired Model for Incremental Learning.  In ICLR, 2018.
  * [22] Gyuhak Kim, Changnan Xiao, Tatsuya Konishi, Zixuan Ke, and Bing Liu.  A theoretical study on solving continual learning.  In Advances in Neural Information Processing Systems, 2022.
  * [23] James Kirkpatrick, Razvan Pascanu, Neil Rabinowitz, Joel Veness, Guillaume Desjardins, Andrei A Rusu, Kieran Milan, John Quan, Tiago Ramalho, Agnieszka Grabska-Barwinska, et al.  Overcoming catastrophic forgetting in neural networks.  Proceedings of the national academy of sciences, 114(13):3521–3526, 2017.
  * [24] Łukasz Korycki and Bartosz Krawczyk.  Class-incremental experience replay for continual learning under concept drift.  arXiv preprint arXiv:2104.11861, 2021.
  * [25] Alex Krizhevsky and Geoffrey Hinton.  Learning multiple layers of features from tiny images.  Technical Report TR-2009, University of Toronto, Toronto., 2009\.
  * [26] Ya Le and Xuan Yang.  Tiny imagenet visual recognition challenge.  CS 231N, 7:7, 2015.
  * [27] Yann LeCun, Corinna Cortes, and Christopher JC Burges.  The mnist database of handwritten digits.  http://yann.lecun.com/exdb/mnist/, 1998.
  * [28] Xilai Li, Yingbo Zhou, Tianfu Wu, Richard Socher, and Caiming Xiong.  Learn to grow: A continual structure learning framework for overcoming catastrophic forgetting.  In ICML, 2019.
  * [29] David Lopez-Paz and Marc’Aurelio Ranzato.  Gradient Episodic Memory for Continual Learning.  In NIPS, pages 6470–6479, 2017.
  * [30] Zheda Mai, Ruiwen Li, Hyunwoo Kim, and Scott Sanner.  Supervised contrastive replay: Revisiting the nearest class mean classifier in online class-incremental continual learning.  In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, pages 3589–3599, 2021.
  * [31] Zheda Mai, Ruiwen Li, Hyunwoo Kim, and Scott Sanner.  Supervised contrastive replay: Revisiting the nearest class mean classifier in online class-incremental continual learning.  In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition (CVPR) Workshops, pages 3589–3599, 2021.
  * [32] Marc Masana, Xialei Liu, Bartlomiej Twardowski, Mikel Menta, Andrew D Bagdanov, and Joost van de Weijer.  Class-incremental learning: survey and performance evaluation on image classification.  arXiv preprint arXiv:2010.15277, 2020.
  * [33] Sudhanshu Mittal, Silvio Galesso, and Thomas Brox.  Essentials for class incremental learning.  In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, pages 3513–3522, 2021.
  * [34] Oleksiy Ostapenko, Mihai Puscas, Tassilo Klein, Patrick Jahnichen, and Moin Nabi.  Learning to remember: A synaptic plasticity driven framework for continual learning.  In CVPR, pages 11321–11329, 2019.
  * [35] Ameya Prabhu, Philip HS Torr, and Puneet K Dokania.  Gdumb: A simple approach that questions our progress in continual learning.  In EECV, pages 524–540, 2020.
  * [36] Jathushan Rajasegaran, Munawar Hayat, Salman Khan, Fahad Shahbaz, and Khan Ling Shao.  Random path selection for incremental learning.  In NeurIPS, 2019.
  * [37] Sylvestre-Alvise Rebuffi, Alexander Kolesnikov, Georg Sperl, and Christoph H Lampert.  icarl: Incremental classifier and representation learning.  In CVPR, pages 2001–2010, 2017.
  * [38] Ari Seff, Alex Beatson, Daniel Suo, and Han Liu.  Continual learning in generative adversarial nets.  arXiv preprint arXiv:1705.08395, 2017.
  * [39] Joan Serrà, Dídac Surís, Marius Miron, and Alexandros Karatzoglou.  Overcoming catastrophic forgetting with hard attention to the task.  In ICML, 2018.
  * [40] Dongsub Shim, Zheda Mai, Jihwan Jeong, Scott Sanner, Hyunwoo Kim, and Jongseong Jang.  Online class-incremental continual learning with adversarial shapley value.  In Proceedings of the AAAI Conference on Artificial Intelligence (AAAI), pages 9630–9638, 2021.
  * [41] Hanul Shin, Jung Kwon Lee, Jaehong Kim, and Jiwon Kim.  Continual learning with deep generative replay.  In NIPS, pages 2994–3003, 2017.
  * [42] Albin Soutif-Cormerais, Marc Masana, Joost Van de Weijer, and Bartlømiej Twardowski.  On the importance of cross-task features for class-incremental learning.  arXiv: 2106.11930, 2021.
  * [43] Johannes von Oswald, Christian Henning, João Sacramento, and Benjamin F Grewe.  Continual learning with hypernetworks.  ICLR, 2020.
  * [44] Liyuan Wang, Xingxing Zhang, Kuo Yang, Longhui Yu, Chongxuan Li, Lanqing Hong, Shifeng Zhang, Zhenguo Li, Yi Zhong, and Jun Zhu.  Memory replay with data compression for continual learning.  ICLR-2022, 2022.
  * [45] Chenshen Wu, Luis Herranz, Xialei Liu, Joost van de Weijer, Bogdan Raducanu, et al.  Memory replay gans: Learning to generate new categories without forgetting.  In NIPS, pages 5962–5972, 2018.
  * [46] Yue Wu, Yinpeng Chen, Lijuan Wang, Yuancheng Ye, Zicheng Liu, Yandong Guo, and Yun Fu.  Large scale incremental learning.  In Proceedings ofthe IEEE/CVF Conference on Computer Vision and Pattern Recognition (CVPR), 2019.
  * [47] Shipeng Yan, Jiale Zhou, Jiangwei Xie, Songyang Zhang, and Xuming He.  An em framework for online incremental learning of semantic segmentation.  In Proceedings of the 29th ACM International Conference on Multimedia, pages 3052–3060, 2021.
  * [48] Haiyan Yin, Ping Li, et al.  Mitigating forgetting in online continual learning with neuron calibration.  Advances in Neural Information Processing Systems, 34, 2021.
  * [49] Jaehong Yoon, Eunho Yang, Jeongtae Lee, and Sung Ju Hwang.  Lifelong Learning with Dynamically Expandable Networks.  In ICLR, 2018.
  * [50] Lu Yu, Bartlomiej Twardowski, Xialei Liu, Luis Herranz, Kai Wang, Yongmei Cheng, Shangling Jui, and Joost van de Weijer.  Semantic drift compensation for class-incremental learning.  In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, pages 6982–6991, 2020.
  * [51] Guanxiong Zeng, Yang Chen, Bo Cui, and Shan Yu.  Continual learning of context-dependent processing in neural networks.  Nature Machine Intelligence, 1(8):364–372, 2019.
  * [52] Junting Zhang, Jie Zhang, Shalin Ghosh, Dawei Li, Serafettin Tasci, Larry Heck, Heming Zhang, and C.-C. Jay Kuo.  Class-incremental learning via deep model consolidation.  In CVPR, 2020.
  * [53] Bowen Zhao, Xi Xiao, Guojun Gan, Bin Zhang, and Shu-Tao Xia.  Maintaining discrimination and fairness in class incremental learning.  In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, pages 13208–13217, 2020.
  * [54] Hanbin Zhao, Hui Wang, Yongjian Fu, Fei Wu, and Xi Li.  Memory efficient class-incremental learning for image classification.  IEEE Transactions on Neural Networks and Learning Systems, 2021\.
  * [55] Fei Zhu, Zhen Cheng, Xu-yao Zhang, and Cheng-lin Liu.  Class-incremental learning via dual augmentation.  Advances in Neural Information Processing Systems, 34, 2021.