# A comprehensive theory of second-order spin photocurrents

Jae-Mo Lihm<EMAIL_ADDRESS>Cheol-Hwan Park<EMAIL_ADDRESS>Center
for Correlated Electron Systems, Institute for Basic Science, Seoul 08826,
Korea  Department of Physics and Astronomy, Seoul National University, Seoul
08826, Korea  Center for Theoretical Physics, Seoul National University, Seoul
08826, Korea

###### Abstract

The spin photocurrents, direct currents induced by light, hold great promise
for introducing new elements to spintronics. However, a general theory for
spin photocurrents in real materials which is applicable to systems with spin-
orbit coupling or noncollinear magnetism is absent. Here, we develop such a
general theory of second-order spin photocurrents. We find that the second-
order spin photocurrents can be classified into Drude, Berry curvature dipole,
shift, injection, and rectification currents, which have different physical
origins and symmetry properties. Surprisingly, our theory predicts a direct
pure spin rectification current in an insulator induced by photons with
energies lower than the material band gap. This phenomenon is absent in the
case of the charge photocurrent. We find that the pure spin current of BiTeI
induced by subgap light is large enough to be observable in experiments.
Moreover, the subgap pure spin photocurrent is highly tunable with the
polarization of light and the flowing direction of the spin photocurrent. This
study lays the groundwork for the study of nonlinear spin photocurrents in
real materials and provides a route to engineer light-controlled spin
currents.

##  I Introduction

Photocurrents, the current response of materials under light irradiation, has
potential in solar cell applications Grinberg  _et al._ (2013) and
optoelectronic devices Bonaccorso  _et al._ (2010). The spin version of the
photocurrent response is gaining much attention as a promising tool in
spintronics Bhat  _et al._ (2005); Ivchenko and Tarasenko (2008); Young  _et
al._ (2013). Spin photocurrents enable an optical generation of spin-polarized
currents, essential building blocks of spintronics. Spin photocurrents have
been intensely investigated, both theoretically Bhat  _et al._ (2005);
Tarasenko and Ivchenko (2005); Young  _et al._ (2013) and experimentally Zhao
_et al._ (2005); Ganichev  _et al._ (2006).

However, a general theory of bulk spin photocurrents is yet to be developed,
despite its importance and the interest in it. Most of the existing studies
either deal with spin-collinear systems, where the spin-up and spin-down
electrons are completely decoupled Young  _et al._ (2013); Fei  _et al._
(2020); Xiao  _et al._ (2021), or consider a specific part of the full
response, such as the Drude current Ganichev  _et al._ (2001); Hamamoto  _et
al._ (2017) or the injection current Bhat  _et al._ (2005); Cui  _et al._
(2007); Fei  _et al._ (2021). The spin shift current in spin-orbit coupled
systems has been investigated in a few studies Bhat  _et al._ (2005); Kim  _et
al._ (2017); Xu  _et al._ (2021); Mu  _et al._ (2021) although we show that
those descriptions for the spin shift current are not complete. In addition,
noncollinear magnetism enables the generation of spin-polarized currents
Železný  _et al._ (2017) and spin Hall currents Zhang  _et al._ (2018) without
spin-orbit coupling. However, current theories cannot describe spin
photocurrents in such systems. This situation is in sharp contrast to the case
of charge photocurrents Sipe and Shkrebtii (2000), where a complete
classification of second-order responses was recently developed Zhang  _et
al._ (2019); de Juan  _et al._ (2020); Holder  _et al._ (2020); Ahn  _et al._
(2020); Gao  _et al._ (2020); Watanabe and Yanase (2021). Thus, a complete,
generally applicable theory of spin photocurrent at the level of its charge
counterpart is highly desirable.

In this paper, we develop a complete theory of second-order spin
photocurrents, which describes the spin and charge photocurrents in a unified
framework. Using perturbation theory in the length gauge, we derive the
expression for the second-order spin conductivity tensor. We show that the
spin photocurrent can be classified into Drude, Berry curvature dipole, shift,
injection, and rectification currents. We find that the spin shift and
rectification currents can be written in terms of the complex-valued spin
shift vector, a gauge-invariant quantity proposed in this work. Interestingly,
the spin rectification current is nonzero in insulators, in sharp contrast
with its charge counterpart, the “intrinsic Fermi surface” current, which is
always zero in insulators de Juan  _et al._ (2020); Gao  _et al._ (2020);
Watanabe and Yanase (2021). Most importantly, finite spin rectification
currents can be generated even by light whose frequency is lower than the band
gap. These subgap pure spin currents are highly tunable with the polarization
of light. We demonstrate our theory and findings by presenting our
calculations on the spin and charge photoconductivity of BiTeI.

##  II Theory

###  II.1 Second-order spin photoconductivity

We study the response of the system to external electric fields using the
Schrödinger equation for the one-particle reduced density matrix Ventura  _et
al._ (2017). We write the responses in terms of the velocity matrix element

$v^{a}_{mn,{\bm{\mathbf{k}}}}=\matrixelement{u_{m\mathbf{k}}}{\partial^{a}H_{0}({\bm{\mathbf{k}}})}{u_{n\mathbf{k}}}/\hbar,$
(1)

occupation factor
$f_{m\mathbf{k}}=1/\\{1+\exp[(\varepsilon_{m\mathbf{k}}-\mu)/k_{\rm B}T]\\}$,
and frequency $\omega_{mn}=(\varepsilon_{m}-\varepsilon_{n})/\hbar$. Here,
$H_{0}({\bm{\mathbf{k}}})$ is the periodic Bloch Hamiltonian,
$\ket{u_{m\mathbf{k}}}$ the periodic part of the Bloch wavefunction of state
$m$ with crystal momentum ${\bm{\mathbf{k}}}$, $\varepsilon_{m\mathbf{k}}$ the
band energy, $T$ the temperature, $\mu$ the chemical potential, and
$\partial^{a}=\partial/\partial k_{a}$. In the following, we omit the
subscript ${\bm{\mathbf{k}}}$ for brevity. We also define
$f_{mn}=f_{m}-f_{n}$. We define the Berry connection
$\xi^{a}_{mn}=i\innerproduct{u_{m}}{\partial^{a}u_{n}}$, which satisfies

$\xi^{a}_{mn}=-iv^{a}_{mn}/\omega_{mn}$ (2)

for non-degenerate states with $m\neq n$.

The spin current is the expectation value of the spin-current operator, whose
matrix element is

$j^{s,a}_{mn}=\frac{1}{2}\matrixelement{u_{m}}{\anticommutator{S^{s}}{v^{a}}}{u_{n}},$
(3)

where $S^{s}$ is the spin operator with $s=x,y,z$. We also let $j^{0,a}=v^{a}$
be the charge-current operator by defining $S^{0}$ to be the identity
operator. For later use, we define

$\Delta j^{s,a}_{mn}=j^{s,a}_{mm}-j^{s,a}_{nn}.$ (4)

(See Appendix A for the discussion on alternate forms of the spin-current
operator.)

Another key quantity for nonlinear spin currents is the spin-velocity
derivative of the velocity operator, which we define as follows:

$d^{s,b;a}_{mn}=j^{s,ab}_{mn}+\sum_{p\neq
m}\frac{j^{s,a}_{mp}v^{b}_{pn}}{\omega_{mp}}+\sum_{p\neq
n}\frac{v^{b}_{mp}j^{s,a}_{pn}}{\omega_{np}}.$ (5)

Here, we also defined

$j^{s,ab}_{mn}=\frac{1}{2}\matrixelement{u_{m}}{\anticommutator{S^{s}}{\partial^{{a}}\partial^{{b}}H_{0}({\bm{\mathbf{k}}})}}{u_{n}}.$
(6)

The name “spin-velocity derivative” reflects the fact that for the charge
current, $d^{0,b;a}_{mn}$ is identical to the generalized derivative of the
velocity operator:

$d^{s=0,b;a}_{mn}=\partial^{{a}}v^{b}_{mn}-i(\xi^{a}_{mm}-\xi^{a}_{nn})v^{b}_{mn}.$
(7)

For spin-collinear systems, the spin-velocity derivatives are equivalent to
the charge counterparts except for a minus sign for the spin-down states. On
the contrary, for spin-noncollinear systems, the spin-velocity derivatives
cannot be simply related to the charge counterparts.

In this work, we focus on the second-order response, which is the lowest order
where a DC current response can occur under AC driving fields. We consider
external fields with frequencies $-\Omega$ and $\Omega+\omega$ and study the
second-order response with frequency $\omega$ in the DC limit
$\omega\rightarrow 0$ with $\Omega$ fixed. We assume a clean system with an
infinite carrier lifetime as done in previous studies of charge currents Sipe
and Shkrebtii (2000); Parker  _et al._ (2019); Zhang  _et al._ (2019); de Juan
_et al._ (2020); Gao  _et al._ (2020); Ahn  _et al._ (2020); Watanabe and
Yanase (2021). Thus, the DC limit should not be taken literally: the
photocurrent should be understood as that from a difference frequency
generation where $\Omega\gg\omega\gg 1/\tau$ de Juan  _et al._ (2020), where
$\tau$ is the characteristic lifetime of the bands.

Watanabe and Yanase Watanabe and Yanase (2021) applied the general results of
Ref. Ventura  _et al._ (2017) to the second-order responses of charge
currents. In this work, we extend this formalism Ventura  _et al._ (2017);
Watanabe and Yanase (2021) to spin currents. The second-order spin and charge
photoconductivity is the sum of five distinct responses: Drude, Berry
curvature dipole (BCD), shift, injection, and rectification
photoconductivities. The formula for each photoconductivity is as follows.

$\displaystyle\sigma^{s,a;bc}_{\mathrm{Drude}}(\Omega)=\frac{q^{3}}{2\hbar^{2}V\Omega^{2}}\sum_{{\bm{\mathbf{k}}},m}j^{s,a}_{mm}\partial^{{b}}\partial^{{c}}f_{m}$
(8)
$\displaystyle\sigma^{s,a;bc}_{\mathrm{BCD}}(\Omega)=\frac{-iq^{3}}{\hbar^{2}V\Omega}\sum_{\begin{subarray}{c}{\bm{\mathbf{k}}},m,n\\\
m\neq
n\end{subarray}}\frac{\imaginary(j^{s,a}_{mn}v^{b}_{nm})}{\omega_{mn}^{2}}\partial^{{c}}f_{m}-(b\leftrightarrow
c)$ (9) $\displaystyle\sigma^{s,a;bc}_{\mathrm{inj.,\tau}}(\Omega)$
$\displaystyle=\tau\eta^{s,a;bc}_{\mathrm{inj.}}(\Omega)$ (10)
$\displaystyle=-\tau\frac{\pi
q^{3}}{\hbar^{2}\Omega^{2}V}\sum_{{\bm{\mathbf{k}}},m,n}f_{mn}\Delta
j^{s,a}_{mn}v^{b}_{mn}v^{c}_{nm}\delta(\Omega+\omega_{mn})$
$\displaystyle\sigma^{s,a;bc}_{\mathrm{shift}}(\Omega)=\frac{i\pi
q^{3}}{2\hbar^{2}V\Omega^{2}}\sum_{{\bm{\mathbf{k}}},m,n}$
$\displaystyle(d^{s,b;a}_{mn}v^{c}_{nm}-d^{s,c;a}_{nm}v^{b}_{mn})$
$\displaystyle\times f_{mn}\delta(\Omega+\omega_{mn})$ (11)
$\displaystyle\sigma^{s,a;bc}_{\mathrm{rect.}}(\Omega)=\frac{q^{3}}{2\hbar^{2}V}\sum_{\begin{subarray}{c}{\bm{\mathbf{k}}},m,n\\\
m\neq
n\end{subarray}}\frac{f_{mn}}{\omega_{mn}^{2}}\Bigg{[}\Big{(}d^{s,b;a}_{mn}v^{c}_{nm}+d^{s,c;a}_{nm}v^{b}_{mn}$
$\displaystyle-\frac{2\Delta
j^{s,a}_{mn}v^{b}_{mn}v^{c}_{nm}}{\omega_{mn}}\Big{)}\mathrm{P}\frac{1}{\Omega+\omega_{mn}}-\frac{\Delta
j^{s,a}_{mn}v^{b}_{mn}v^{c}_{nm}}{(\Omega+\omega_{mn})^{2}}\Bigg{]}$ (12)

Here, $V$ is the volume of the system, and $q$ the charge of an electron. The
detailed derivation of these equations is given in Appendix B. For the
injection current, we introduced a phenomenological relaxation time $\tau$
Passos  _et al._ (2018). This treatment is needed to avoid divergence in the
DC limit and is justified by calculations based on the Floquet formalism de
Juan  _et al._ (2017).

The Drude [Eq. (8)] and BCD [Eq. (9)] terms contain derivatives of the
occupation factor. Therefore, the corresponding currents are zero in
insulators. While the shift and injection currents originate from the
absorptive (resonant) response, the rectification current originates from the
reactive (non-resonant) response. The injection and shift currents are nonzero
only when occupied and unoccupied states are resonantly coupled by light:
$\Omega+\omega_{mn}=0$. In contrast, the rectification current does not
require an energy-conserving transition.

The charge rectification current was discovered only recently de Juan  _et
al._ (2020); Gao  _et al._ (2020); Watanabe and Yanase (2021). It was termed
the “intrinsic Fermi surface” contribution because one can convert Eq. (II.1)
for the charge current case into a Fermi-surface integral using Eq. (7):

$\displaystyle\sigma^{i=0,{a;bc}}_{\mathrm{rect.}}=\frac{q^{3}}{2\hbar^{2}V}\sum_{{\bm{\mathbf{k}}},m,n}f_{mn}\partial^{{a}}\left(\frac{v^{b}_{mn}v^{c}_{nm}}{\omega_{mn}^{2}}\mathrm{P}\frac{1}{\Omega+\omega_{mn}}\right)$
$\displaystyle=-\frac{q^{3}}{2\hbar^{2}V}\sum_{{\bm{\mathbf{k}}},m,n}(\partial^{{a}}f_{mn})\frac{v^{b}_{mn}v^{c}_{nm}}{\omega_{mn}^{2}}\mathrm{P}\frac{1}{\Omega+\omega_{mn}}.$
(13)

Hence, the charge rectification current is zero in insulators.

However, we find that the spin rectification current is nonzero even in
insulators. One of the reasons for the contrast with the charge case is that
the spin-velocity derivative $d^{s,b;a}_{mn}$ with $s\neq 0$ is not an actual
derivative with respect to the crystal momentum. In other words, the spin-
velocity derivative, like the spin Berry curvature Guo  _et al._ (2008), is
not a geometric quantity. Thus, in the spin current case, one cannot rewrite
$\sigma^{s,{a;bc}}_{\mathrm{rect.}}$ as a Fermi-surface integral.

Since the rectification current does not require a resonant coupling of
occupied and unoccupied states, it leads to a subgap spin photocurrent: a DC
spin photocurrent under irradiation of light with a frequency less than the
band gap of an insulating system. Such a subgap response is unique to the spin
photocurrent, while absent in the charge counterpart.

We note that the previously reported equations for the spin shift current Kim
_et al._ (2017); Xu  _et al._ (2021); Mu  _et al._ (2021) are not equivalent
to our results [Eqs. (II.1, II.2)]. Reference Kim  _et al._ (2017), which uses
the velocity gauge, considers only the diamagnetic two-photon output vertex
and ignores the contribution of a diagram with three one-photon vertices. The
photoconductivity formula of Refs. Xu  _et al._ (2021); Mu  _et al._ (2021),
which is also based on the velocity gauge, does not take the two-photon vertex
into account. To obtain the correct result using the velocity gauge formalism,
one needs to consider all diagrams consisted of one-, two-, and three-photon
vertices and apply related sum rules Passos  _et al._ (2018); Parker  _et al._
(2019); Gao  _et al._ (2020).

###  II.2 Writing photoconductiviy in terms of the complex shift vector

To understand the physical mechanism behind the spin shift and rectification
currents, we now show that these currents can be written in terms of the “spin
shift vector” which is a gauge-invariant quantity we introduce in this work.

For incident light with polarization vector ${\bm{\mathbf{{\mathcal{E}}}}}$,
the second-order DC current response reads

$\sigma^{s,a;{\mathcal{E}}}(\Omega)=\sum_{b,c}\sigma^{s,a;bc}(\Omega){\mathcal{E}}^{b}{\mathcal{E}}^{*c}.$
(14)

We also define

$\mathcal{O}^{\mathcal{E}}_{mn}=\sum_{a}\mathcal{O}^{a}_{mn}{\mathcal{E}}^{a}$
(15)

for any vector or tensor operator ${\bm{\mathbf{\mathcal{O}}}}$. For a
Hermitian matrix $\mathcal{O}$, the following holds:

$\left(\mathcal{O}^{{\mathcal{E}}^{*}}_{nm}\right)^{*}=\sum_{a}\left(\mathcal{O}^{a}_{nm}{\mathcal{E}}^{*a}\right)^{*}=\sum_{a}\mathcal{O}^{a}_{mn}{\mathcal{E}}^{a}=\mathcal{O}^{{\mathcal{E}}}_{mn}.$
(16)

Following Ref. Ahn  _et al._ (2020), we define the complex shift vector
between two different bands $m$ and $n$ for light polarization
${\bm{\mathbf{{\mathcal{E}}}}}$ as

${\widetilde{R}}^{{\mathcal{E}};a}_{mn}=i\frac{\mathcal{D}^{a}\xi^{\mathcal{E}}_{mn}}{\xi^{\mathcal{E}}_{mn}}=\xi^{a}_{mm}-\xi^{a}_{nn}+i\partial^{a}\log\xi^{\mathcal{E}}_{mn}$
(17)

Note that we consider a fixed light polarization because the shift vector does
not transform like a vector with respect to the rotation of the light
polarization. By using Eq. (2), the complex shift vector can be written in
terms of the velocity matrix and its generalized derivative:

${\widetilde{R}}^{{\mathcal{E}};a}_{mn}=i\left(\frac{d^{0,{\mathcal{E}};a}_{mn}}{v^{\mathcal{E}}_{mn}}-\frac{\Delta
v^{a}_{mn}}{\omega_{mn}}\right).$ (18)

The standard shift vector Sipe and Shkrebtii (2000) is the real-part of the
complex shift vector:

$R^{{\mathcal{E}};a}_{mn}=\xi^{a}_{mm}-\xi^{a}_{nn}-\partial^{a}\arg\xi^{\mathcal{E}}_{mn}=\real{\widetilde{R}}^{{\mathcal{E}};a}_{mn},$
(19)

where ‘$\arg$’ denotes the argument of a complex number. This real-valued
shift vector describes the change in the position of the electron in the
transition from band $n$ to $m$.

Now, we define the spin shift vector as the spin generalization of Eq. (18):

${\widetilde{R}}^{s,{\mathcal{E}};a}_{mn}=i\left(\frac{d^{s,{\mathcal{E}};a}_{mn}}{v^{\mathcal{E}}_{mn}}-\frac{\Delta
j^{s,a}_{mn}}{\omega_{mn}}\right).$ (20)

The spin shift vector is invariant to the gauge transformation
$\ket{u_{m\mathbf{k}}}\rightarrow\ket{u_{m\mathbf{k}}}e^{i\phi_{m\mathbf{k}}}$
and thus represents a physical quantity. By an analogy to the charge case, the
real part of the spin shift vector can be interpreted as the shift in the spin
position in the transition from band $n$ to $m$. However, in contrast to the
charge case, it is not possible to write
${\widetilde{R}}^{s,{\mathcal{E}};a}_{mn}$ directly in terms of the matrix
elements between bands $m$ and $n$ and their derivatives due to the
nongeometric nature of the spin-current operator.

To rewrite the shift photoconductivity [Eq. (II.1)] and the rectification
conductivity [Eq. (II.1)] in terms of the spin shift vector, we use the
following property:

$\displaystyle\left({\widetilde{R}}^{s,{\mathcal{E}};a}_{mn}\right)^{*}=$
$\displaystyle-i\left(\frac{d^{s,{\mathcal{E}};a}_{mn}}{v^{{\mathcal{E}}}_{mn}}-\frac{\Delta
j^{s,a}_{mn}}{\omega_{mn}}\right)^{*}$ $\displaystyle=$
$\displaystyle-i\left(\frac{d^{s,{\mathcal{E}}^{*};a}_{nm}}{v^{{\mathcal{E}}^{*}}_{nm}}-\frac{\Delta
j^{s,a}_{nm}}{\omega_{nm}}\right)$ $\displaystyle=$
$\displaystyle-{\widetilde{R}}^{s,{\mathcal{E}}^{*};a}_{nm}.$ (21)

Then, by substituting Eq. (II.1) into Eq. (14), we can write the spin shift
current for light polarization ${\mathcal{E}}$ in terms of the spin shift
vector as

$\displaystyle\sigma^{s,a;{\mathcal{E}}}_{\mathrm{shift}}(\Omega)=$
$\displaystyle\frac{i\pi
q^{3}}{2\hbar^{2}V\Omega^{2}}\sum_{{\bm{\mathbf{k}}},m,n}(d^{s,{\mathcal{E}};a}_{mn}v^{{\mathcal{E}}^{*}}_{nm}-d^{s,{\mathcal{E}}^{*};a}_{nm}v^{\mathcal{E}}_{mn})f_{mn}\delta(\Omega+\omega_{mn})$
$\displaystyle=$ $\displaystyle\frac{i\pi
q^{3}}{2\hbar^{2}V\Omega^{2}}\sum_{{\bm{\mathbf{k}}},m,n}(-i{\widetilde{R}}^{s,{\mathcal{E}};a}_{mn}v^{\mathcal{E}}_{mn}v^{{\mathcal{E}}^{*}}_{nm}+\Delta
j^{s,a}_{mn}\frac{v^{\mathcal{E}}_{mn}v^{{\mathcal{E}}^{*}}_{nm}}{\omega_{mn}}+i{\widetilde{R}}^{s,{\mathcal{E}}^{*};a}_{nm}v^{{\mathcal{E}}^{*}}_{nm}v^{\mathcal{E}}_{mn}-\Delta
j^{s,a}_{nm}\frac{v^{{\mathcal{E}}^{*}}_{nm}v^{\mathcal{E}}_{mn}}{\omega_{nm}})f_{mn}\delta(\Omega+\omega_{mn})$
$\displaystyle=$ $\displaystyle\frac{i\pi
q^{3}}{2\hbar^{2}V\Omega^{2}}\sum_{{\bm{\mathbf{k}}},m,n}(-i{\widetilde{R}}^{s,{\mathcal{E}};a}_{mn}+i{\widetilde{R}}^{s,{\mathcal{E}}^{*};a}_{nm})f_{mn}\omega_{mn}^{2}\xi^{\mathcal{E}}_{mn}\xi^{{\mathcal{E}}^{*}}_{nm}\delta(\Omega+\omega_{mn})$
$\displaystyle=$ $\displaystyle\frac{\pi
q^{3}}{\hbar^{2}V}\sum_{{\bm{\mathbf{k}}},m,n}\real\left({\widetilde{R}}^{s,{\mathcal{E}};a}_{mn}\right)f_{mn}\absolutevalue{\xi^{\mathcal{E}}_{mn}}^{2}\delta(\Omega+\omega_{mn}).$
(22)

In the last equality, we used Eq. (II.2). For the charge case $s=0$, Eq.
(II.2) agrees with the known formula for the shift current under linear and
circular polarization Sipe and Shkrebtii (2000); Ahn  _et al._ (2020);
Watanabe and Yanase (2021). (Note that according to our derivation, the shift
current due to circularly polarized light corresponds to the “gyration
current” term in Ref. Watanabe and Yanase (2021).) The novel finding of this
work is that the spin shift current can be written in the same form by
defining the spin shift vector as Eq. (20).

The last line of Eq. (II.2) has the form of the Fermi golden rule. Hence, the
spin shift current can be interpreted as the current due to the spin density
dipole generated by the transition from state $n$ to $m$. This interpretation
parallels the case of charge shift current von Baltz and Kraut (1981); Sipe
and Shkrebtii (2000). However, in contrast to the charge case [Eq. (19)], it
is not possible to write the spin shift vector
${\widetilde{R}}^{s,{\mathcal{E}};a}_{mn}$ only in terms of the low-energy
wavefunctions due to the nongeometric nature of the spin-current operator.

For the spin rectification current, we find

$\displaystyle\sigma^{s,a;{\mathcal{E}}}_{\mathrm{rect.}}(\Omega)=$
$\displaystyle\frac{q^{3}}{2\hbar^{2}V}\sum_{\begin{subarray}{c}{\bm{\mathbf{k}}},m,n\\\
m\neq
n\end{subarray}}\frac{f_{mn}}{\omega_{mn}^{2}}\Bigg{[}\Big{(}d^{s,{\mathcal{E}};a}_{mn}v^{{\mathcal{E}}^{*}}_{nm}+d^{s,{\mathcal{E}}^{*};a}_{nm}v^{\mathcal{E}}_{mn}-\frac{2\Delta
j^{s,a}_{mn}v^{\mathcal{E}}_{mn}v^{{\mathcal{E}}^{*}}_{nm}}{\omega_{mn}}\Big{)}\mathrm{P}\frac{1}{\Omega+\omega_{mn}}-\frac{\Delta
j^{s,a}_{mn}v^{\mathcal{E}}_{mn}v^{{\mathcal{E}}^{*}}_{nm}}{(\Omega+\omega_{mn})^{2}}\Bigg{]}$
$\displaystyle=$
$\displaystyle\frac{q^{3}}{2\hbar^{2}V}\sum_{\begin{subarray}{c}{\bm{\mathbf{k}}},m,n\\\
m\neq
n\end{subarray}}\frac{f_{mn}}{\omega_{mn}^{2}}\Bigg{[}\left(-i{\widetilde{R}}^{s,{\mathcal{E}};a}_{mn}v^{\mathcal{E}}_{mn}v^{{\mathcal{E}}^{*}}_{nm}-i{\widetilde{R}}^{s,{\mathcal{E}}^{*};a}_{nm}v^{{\mathcal{E}}^{*}}_{nm}v^{\mathcal{E}}_{mn}\right)\mathrm{P}\frac{1}{\Omega+\omega_{mn}}-\frac{\Delta
j^{s,a}_{mn}v^{\mathcal{E}}_{mn}v^{{\mathcal{E}}^{*}}_{nm}}{(\Omega+\omega_{mn})^{2}}\Bigg{]}$
$\displaystyle=$
$\displaystyle\frac{q^{3}}{\hbar^{2}V}\sum_{\begin{subarray}{c}{\bm{\mathbf{k}}},m,n\\\
m\neq
n\end{subarray}}f_{mn}\Bigg{[}\imaginary\left({\widetilde{R}}^{s,{\mathcal{E}};a}_{mn}\right)\mathrm{P}\frac{1}{\Omega+\omega_{mn}}-\frac{\Delta
j^{s,a}_{mn}}{2(\Omega+\omega_{mn})^{2}}\Bigg{]}\absolutevalue{\xi^{\mathcal{E}}_{mn}}^{2}.$
(23)

In the last equality, we used Eq. (II.2).

Interestingly, the shift and rectification currents are proportional to the
real and imaginary parts of the shift vector, respectively. While the real
part of the shift vector has been understood as the positional shift of the
electron wave packet von Baltz and Kraut (1981); Sipe and Shkrebtii (2000), we
find here that the imaginary part plays an important role in the charge and
spin rectification currents. The physical meaning of the imaginary part
remains a subject of future study.

###  II.3 Symmetry of photoconductivity

| Current type | Jahn symbols
---|---|---
Linear pol. | Curcular pol.
Charge $\sigma^{0,a;bc}$ | Drude | a[V3] | -
Berry curvature dipole | - | V{V2}
Shift | V[V2] | aV{V2}
Injection, Rectification | aV[V2] | V{V2}
Spin $\sigma^{s,a;bc}$ | Drude | eVV[V2] | -
Berry curvature dipole | - | aeVV{V2}
Shift | aeVV[V2] | eVV{V2}
Injection, Rectification | eVV[V2] | aeVV{V2}
Table 1:  Jahn symbols for the second-order charge and spin
photoconductivities in the clean limit. Hyphens indicate that the
corresponding currents are always zero.

Next, we study the symmetry properties of the photocurrents. To represent the
symmetry properties, we use the Jahn symbols Jahn (1949); Gallego  _et al._
(2019). Let us briefly explain the Jahn symbols. Each ‘V’ term corresponds to
the vector indices. The number of ‘V’ corresponds to the tensor rank. V’s
inside [ ] and { } denote symmetric and anti-symmetric indices, respectively.
Symbol ‘e’ and ‘a’ denotes that the tensor is axial and odd under time
reversal, respectively.

The real symmetric (imaginary anti-symmetric) component of the
photoconductivity tensor, $\real(\sigma^{s,a;bc}+\sigma^{s,a;cb})$
($\imaginary(\sigma^{s,a;bc}-\sigma^{s,a;cb})$) corresponds to the generation
of a current under linearly (circularly) polarized light. In Appendix C, we
derive the symmetry transformation properties of the five charge
photoconductivity tensors using the symmetry of the velocity matrix elements.

One can easily derive the symmetry of the spin photoconductivity from the
symmetry of the corresponding charge photoconductivity. The second-order spin
photoconductivity $\sigma^{s,a;bc}$ is a rank-4 tensor due to the additional
vector describing spin polarization. Equation (3) shows that the spin velocity
transforms under symmetry operations like the product of the charge velocity
and the spin polarization vector. The same relation holds for the spin and
charge conductivity. The spin polarization vector is a time-reversal odd,
axial vector. Thus, the Jahn symbol for the spin photoconductivity can be
obtained by multiplying ‘aeV’ to the Jahn symbol of the corresponding charge
photoconductivity. The only exception is the Drude current. As explained in
Appendix C.1, the charge Drude current has an additional symmetry for swapping
the current and field direction indices [Eq. (83)]. An analogous expression
does not exist for the spin Drude current because the spin velocity
$j^{s,a}_{mm}$ for $s=x,y,z$ is not a derivative of some quantity with respect
to the crystal momenta.

In Table 1, we summarize the symmetry of the charge and spin photoconductivity
tensors. The nonzero, independent coefficients for the photoconductivity
tensors of a given material can be easily identified using the MTensor program
Gallego  _et al._ (2019).

Let us focus on the role of the spatial inversion ${\hat{\mathcal{P}}}$, time
reversal ${\hat{\mathcal{T}}}$, and space-time inversion
${\hat{\mathcal{P}}\hat{\mathcal{T}}}$. The properties of the second-order
charge conductivity in inversion- or time-reversal-symmetric systems were
recently studied Ahn  _et al._ (2020); Watanabe and Yanase (2021). Here, we
extend the analysis to the spin conductivity.

Spatial inversion ${\hat{\mathcal{P}}}$ flips the sign of odd-rank vectors and
gives an additional minus sign to axial tensors. Since the charge
photoconductivity is an odd-rank polar (non-axial) tensor, it obtains a minus
sign under ${\hat{\mathcal{P}}}$. The spin photoconductivity is an even-rank
axial tensor, so it also obtains a minus sign under ${\hat{\mathcal{P}}}$.
Therefore, in an inversion-symmetric system, both the spin and charge second-
order photoconductivity is always zero.

Time-reversal operation ${\hat{\mathcal{T}}}$ flips the sign of time-reversal-
odd tensors, which have ‘a’ in their Jahn symbols. Therefore, among the
second-order photoconductivities, only those without ‘a’ in their Jahn symbols
in Table 1 are allowed in a ${\hat{\mathcal{T}}}$-symmetric system. Therefore,
in a ${\hat{\mathcal{T}}}$-symmetric system, the spin (charge) shift current
is induced by circularly (linearly) polarized light, and the spin (charge)
injection and rectification current is induced by linearly (circularly)
polarized light. Also, the spin Drude (charge BCD) current is induced by
linearly (circularly) polarized light.

For a ${\hat{\mathcal{P}}\hat{\mathcal{T}}}$-symmetric system, the relation is
the opposite: only the photocurrents with ‘a’ in their Jahn symbols are
allowed. The reason is as follows. All second-order spin and charge
photoconductivities are all odd under ${\hat{\mathcal{P}}}$. For the
${\hat{\mathcal{P}}\hat{\mathcal{T}}}$ symmetry operation, the $(-1)$ factor
from ${\hat{\mathcal{P}}}$ is multiplied from the sign factor coming from
${\hat{\mathcal{T}}}$. Thus, for both the spin and charge cases, if the
photoconductivity is even (odd) under ${\hat{\mathcal{T}}}$, it is odd (even)
under ${\hat{\mathcal{P}}\hat{\mathcal{T}}}$. In a
${\hat{\mathcal{P}}\hat{\mathcal{T}}}$-symmetric system, only tensors that are
even under ${\hat{\mathcal{P}}\hat{\mathcal{T}}}$ can be nonzero. Therefore,
among the second-order photoconductivities, only those with ‘a’ in their Jahn
symbols (i.e., those which were not allowed in ${\hat{\mathcal{T}}}$-symmetric
system systems) are allowed in a
${\hat{\mathcal{P}}\hat{\mathcal{T}}}$-symmetric system.

The injection and rectification currents have the same symmetry properties. To
distinguish these contributions, one can use a subgap light, which induces
only the rectification current. Another way is to measure the spin current for
samples with different carrier lifetimes (e.g. due to different impurity
concentrations) and use the fact that only the injection current is
proportional to the carrier lifetime.

##  III Results

Figure 1:  (a) Crystal structure of BiTeI. (b) Brillouin zone of BiTeI. The
crystal structure is drawn using the software VESTA Momma and Izumi (2011).
Figure 2:  Electronic energy band structure of BiTeI calculated using density
functional theory (DFT) and Wannier interpolation. The green horizontal dotted
lines indicate the center of the band gap.

Now, we apply our theory to the charge and spin photocurrent response of bulk
BiTeI. BiTeI is a non-magnetic polar direct-gap semiconductor with a giant
Rashba-type spin splitting Ishizaka  _et al._ (2011). Here, we calculate the
charge and spin photoconductivity of BiTeI using a Wannier-function-based ab
initio tight-binding model Marzari and Vanderbilt (1997). Figure 1 shows the
crystal structure and Brillouin zone of BiTeI, and Fig. 2 shows the band
structure of BiTeI.

Figure 3:  (a) Schematic illustration of light-induced spin currents in BiTeI.
(b-d) Charge and (e-g) spin $S^{y}$ photoconductivity of BiTeI. The vertical
dashed lines indicate the band gap of 0.35 eV. The conductivity for linearly
and circularly polarized light corresponds to the real symmetric and imaginary
anti-symmetric components of the conductivity tensor, respectively.

In Fig. 3, we show the charge and spin $S^{y}$ photoconductivity. Spin $S^{x}$
and $S^{z}$ photoconductivities are shown in Fig. 5. Since BiTeI is an
insulator, the charge rectification current is zero in the entire frequency
range. In contrast, the spin rectification current is nonzero, even for subgap
frequencies. Therefore, our theory predicts pure spin currents in the subgap
frequency regime. We note that the spin rectification current diverges when
$\Omega$ approaches the band gap in the limit of zero smearing. In reality,
this divergence will be regularized by the finite broadening of the bands
arising from the finite lifetime of the electronic states.

The $\sigma^{y,x;yy}_{\mathrm{rect.}}$ component of the subgap spin
photoconductivity, which describes the spin $S^{y}$ current along the $x$
direction with irradiation of light linearly polarized along the $y$
direction, has a magnitude around $20~{}\mathrm{\mu A/V^{2}}\hbar/e$ for an
infrared light with frequency 0.30 eV. The magnitude of this subgap spin
response is larger than the calculated spin shift currents of collinear
antiferromagnets BiFeO3 and hematite Young  _et al._ (2013). This pure spin
current could also be detected by using spin-to-charge conversion methods such
as the inverse spin Hall effect Saitoh  _et al._ (2006).

Figure 4:  Spin polarization of the subgap spin rectification conductivity
$\sigma^{s,a}$ [Eq. (24)] of BiTiI. The arrows and colors indicate the in-
plane and out-of-plane spin polarization of the current, respectively. Light
frequency is $\Omega=0.30$ eV, below the calculated band gap of 0.35 eV.

Remarkably, the polarization of the spin current can be tuned by the
polarization of the light. Figure 4 shows the spin polarization of the subgap
spin conductivity as a function of light polarization angles:

$\sigma^{s,a}(\theta,\phi)=\sum_{b,\,c}\sigma^{s,a;bc}\,\hat{n}^{b}\,\hat{n}^{c}\,,$
(24)

where $\hat{n}=(\sin\theta\sin\phi,\,\sin\theta\cos\phi,\,\cos\theta)$ denotes
the direction of the light polarization. The energy of the incident photon is
$\Omega=0.30$ eV, below the calculated band gap of 0.35 eV. Both the direction
and magnitude of the spin polarization can be highly tuned by changing the
light polarization. Also, the spin polarization depends significantly on the
direction of the spin current. These remarkable tunabilities are not only
scientifically important by themselves but may also open doors to novel
spintronics applications.

##  IV Discussion

For a material to produce a large spin photocurrent driven by spin-orbit
coupling, it should consist of atoms with high atomic numbers, have a large
structural asymmetry Tan  _et al._ (2016), and have a small band gap Cook  _et
al._ (2017). BiTeI, which is well known for the large Rashba effect and has a
band gap of around 0.3 eV, is one of the materials satisfying all these
criteria. Monolayer SnS and SnSe, which contain the heavy Sn atoms and host
large charge photocurrents Rangel  _et al._ (2017), are also good candidates.
One could search for other materials that can host large spin photocurrents
based on these basic principles.

In this work, we considered the clean limit and found that only the spin
injection current depends on the lifetime. It has also been numerically tested
that the relaxation time indeed has little effect on the spin and charge
photocurrents Zhang  _et al._ (2019); Xu  _et al._ (2021). The effect of
various scattering processes Ivchenko  _et al._ (1990); Ganichev  _et al._
(2002, 2006); Budkin and Tarasenko (2020) beyond the simple phenomenological
relaxation time approximation on spin photocurrents may be a subject for a
future study.

Recently, Kaplan et al. Kaplan  _et al._ (2020) proposed that a subgap charge
photocurrent exists in systems without time-reversal symmetry. Within our
theory, the subgap charge photocurrent is always zero, with or without the
time-reversal symmetry. We first note that since BiTeI has a time-reversal
symmetry, the subgap charge photocurrent is zero in both theories. The
prediction of subgap photocurrents in a time-reversal symmetric system is
unique to the spin-current response. Also, in this work, we first take the
clean limit of infinite quasiparticle lifetime and then take the DC limit
$\omega\rightarrow 0$ as in the previous studies Sipe and Shkrebtii (2000);
Parker  _et al._ (2019); Zhang  _et al._ (2019); de Juan  _et al._ (2020); Gao
_et al._ (2020); Ahn  _et al._ (2020); Watanabe and Yanase (2021), while Ref.
Kaplan  _et al._ (2020) takes the clean limit after the DC limit. Thus, the
two theories describe different physical situations. Concretely, our results
should be understood as describing the difference frequency generation
scenario, where $\omega\gg 1/\tau$ de Juan  _et al._ (2020). Also, while the
subgap charge photocurrent of Ref. Kaplan  _et al._ (2020) depends on the
detail of how the clean limit is taken, the subgap spin rectification
photocurrent that we report is independent of such detail. The relation
between the subgap spin rectification current and the subgap charge
photocurrent of Ref. Kaplan  _et al._ (2020) is a subject of future study.

##  V Conclusion

In conclusion, we developed a complete theory of second-order spin
photocurrents, which is generally applicable to systems with spin-orbit
coupling or noncollinear magnetization. The spin shift vector, which is a
gauge-invariant quantity, governs the spin shift and rectification responses.
We found that subgap light can induce a DC spin-current response in a
nonmagnetic insulator, which does not occur in the charge current case. This
phenomenon enables an optical generation of pure spin currents, whose spin
polarization is highly tunable by the polarization of light or the flow
direction of the spin current. We applied our theory to BiTeI and found that
the calculated subgap spin photoconductivity is large enough to be detectable
in experiments. Our study paves the way for theoretical and experimental
studies of the nonlinear optical generation of spin currents.

Note added in proof — Calculation of the spin shift [Eq. (II.1)] and
rectification current [Eq. (II.1)], as well as the spin shift vector [Eq.
(20)], requires a sum over an infinite number of bands. In this work, we used
the diagonal tight-binding approximation, which involves an artificial
truncation of the bands. Reference Lihm and Park (2021) presents an efficient
method to calculate these quantities from first principles without any band-
truncation error.

###### Acknowledgements.

This work was supported by the Creative-Pioneering Research Program through
Seoul National University, Korean NRF No-2020R1A2C1014760, and the Institute
for Basic Science (No. IBSR009-D1). Computational resources have been provided
by KISTI (KSC-2020-INO-0078).

##  Appendix A Definition of the spin current

In presence of spin-orbit coupling, spin is not a conserved quantity. Thus,
the conventional definition of the spin-current operator [Eq. (3)], which is
used in most studies including ours, does not describe a conserved current.
Therefore, it is actually still debated whether this conventional, most widely
used definition of the spin-current operator should be modified or not.

A modified definition of the spin-current operator as
$j^{s,a}=d(r^{a}S^{s})/dt$ has been proposed for an inversion-symmetric system
in Ref. Shi  _et al._ (2006). However, this definition gives a conserved
current only if the “spin generation in the bulk is absent.” In the case of
light irradiation, which is a prerequisite for photocurrent response, the
light-induced spin polarization is always present Xu  _et al._ (2021)
Moreover, the spin generation in the bulk through spin-orbit coupling is not
forbidden in general when the bulk is inversion asymmetric, which is another
prerequisite for having a second-order photocurrent response. Thus, even this
modified definition of spin current does not give a conserved current in the
study of photocurrents.

Furthermore, there are some theoretical works claiming that the conventional
definition should be used as is Sun  _et al._ (2008) because the non-
conservation of spin current is a physical phenomenon with possible
experimental outcomes. Also, it has been estimated numerically that the
difference between the conventional and modified spin currents will be on the
order of 10 to 20%, indicating that the computational results will not be
qualitatively affected by the choice of the spin-current operator Xu  _et al._
(2021). Most importantly, the conventional definition has been tested against
experiments in the context of spin Hall effects and gave good agreement on the
spin Hall angle Qiao  _et al._ (2018); Wang  _et al._ (2020). Therefore, we
used the conventional spin-current operator in our work.

We also note that to our knowledge all previous works on the spin
photocurrents, including Refs. Bhat  _et al._ (2005); Fei  _et al._ (2021); Xu
_et al._ (2021), also used the conventional spin-current operator.

One interesting property of the conventional spin-current operator is that a
nonzero spin current exists in equilibrium Rashba (2003). We note in passing
that this equilibrium spin current is an intrinsic property of the spin
current Sun  _et al._ (2008). As we study only the spin currents that occur in
response to the external light, we did not consider the equilibrium spin
current.

##  Appendix B Derivation of spin photoconductivity

In this section, we derive the spin photocurrent conductivity. We assume non-
degenerate bands. (For ${\hat{\mathcal{P}}}{\hat{\mathcal{T}}}$-symmetric
cases, the recipe in Appendix B of Ref. Watanabe and Yanase (2021) can be
used.)

###  B.1 Perturbative expansion of density matrix

Under the independent particle approximation, the electron Hamiltonian reads

$\hat{H}_{0}=\sum_{{\bm{\mathbf{k}}},m}\varepsilon_{m\mathbf{k}}\hat{c}^{\dagger}_{m\mathbf{k}}\hat{c}_{m\mathbf{k}}.$
(25)

Here, $\varepsilon_{m\mathbf{k}}$ is the band energy of state $m$ with crystal
momentum ${\bm{\mathbf{k}}}$, and $\hat{c}^{\dagger}_{m\mathbf{k}}$ and
$\hat{c}_{m\mathbf{k}}$ are the electron creation and annihilation operators,
respectively. Using the Bloch theorem, one can write an eigenvalue equation
for the band energy and the periodic part of the Bloch wavefunction:

$H_{0}({\bm{\mathbf{k}}})\ket{u_{m\mathbf{k}}}=\varepsilon_{m\mathbf{k}}\ket{u_{m\mathbf{k}}}.$
(26)

In the following, we omit the subscript ${\bm{\mathbf{k}}}$ for brevity.

In the length gauge and the dipole approximation, the perturbation due to the
electric field is given by

$\hat{V}(t)=-q\hat{r}^{a}E^{a}(t),$ (27)

where $q$ is the charge of an electron, and $E^{a}$ and $r^{a}$ are the
electric field and position operator along the direction $a$, respectively.
The sum over repeated superscripts is implied. Here, $\hat{r}^{a}$ is the
position operator, whose matrix element is Blount (1962)

$r^{a}_{mn}=i\delta_{mn}\partial^{a}+\xi^{a}_{mn}$ (28)

where we defined $\partial^{a}=\partial/\partial k_{a}$ and the Berry
connection $\xi^{a}_{mn}=i\innerproduct{u_{m}}{\partial^{a}u_{n}}$.

The Schrödinger equation for the one-particle reduced density matrix
$\rho_{mn}$ is Ventura  _et al._ (2017)

$i\hbar\frac{d\rho_{mn}(t)}{dt}=\hbar\omega_{mn}\rho_{mn}(t)-qE^{a}(t)[r^{a},\rho(t)]_{mn}.$
(29)

Defining the Fourier transformation as

$\rho_{ab}(t)=\int\frac{d\omega}{2\pi}e^{-i\omega t}\rho_{ab}(\omega),$ (30)

the Schrödinger equation in the frequency domain becomes

$\hbar(\omega-\omega_{mn})\rho_{mn}(\omega)=-q\int\frac{d\Omega}{2\pi}E^{a}(\Omega)\commutator{r^{a}}{\rho(\omega-\Omega)}_{mn}.$
(31)

Regarding the electric field as a perturbation, one can expand $\rho$ as

$\rho(\omega)=\sum_{n}\rho^{(n)}(\omega),$ (32)

where $\rho^{(n)}$ is of order $\mathcal{O}(E^{n}).$ The zeroth-order term is
the Fermi-Dirac distribution:

$\rho^{(0)}_{ab}(\omega)=2\pi\delta(\omega)f_{a}\delta_{ab}.$ (33)

Using the perturbative expansion of the reduced density matrix, the $n$-th
order charge and spin current responses can be calculated as

$\displaystyle
J^{s,a}_{(n)}(\omega)=\frac{1}{V}\sum_{{\bm{\mathbf{k}}},m,n}qj^{s,a}_{mn}\rho^{(n)}_{nm}(\omega).$
(34)

Here, $V$ is the volume of the system.

In this work, we focus on the second-order response, which is the lowest order
where a DC current response can occur under AC driving fields. We define the
second-order conductivity tensor
$\sigma^{s,a;bc}{(\omega;\omega_{1},\omega_{2})}$ as

$\displaystyle
J^{s,a}_{(2)}(\omega)=\int\frac{d\omega_{1}d\omega_{2}}{(2\pi)^{2}}$
$\displaystyle\sigma^{s,a;bc}{(\omega;\omega_{1},\omega_{2})}E^{b}(\omega_{1})E^{c}(\omega_{2})$
$\displaystyle\times 2\pi\delta(\omega-\omega_{1}-\omega_{2}).$ (35)

To investigate the DC photocurrent, we set

$\omega_{1}=-\Omega,\quad\omega_{2}=\Omega+\omega$ (36)

and take the DC limit $\omega\rightarrow 0$. In other words, we study
$\sigma^{s,a;bc}(0;-\Omega,\Omega)$. We assume a clean system with an infinite
carrier lifetime as done in previous studies on charge-current responses Sipe
and Shkrebtii (2000); Parker  _et al._ (2019); Zhang  _et al._ (2019); de Juan
_et al._ (2020); Gao  _et al._ (2020); Ahn  _et al._ (2020); Watanabe and
Yanase (2021).

Watanabe and Yanase Watanabe and Yanase (2021) showed that the matrix elements
of the second-order reduced density operator can be divided into terms
originating from the intraband (i) and interband (e) position operators as

$\rho_{mn}^{(2)}=\rho_{mn}^{\rm(ii)}+\rho_{mn}^{\rm(ei)}+\rho_{mn}^{\rm(ie)}+\rho_{mn}^{\rm(ee)}.$
(37)

Following Ref. Watanabe and Yanase (2021), we define the contribution of each
term of Eq. (37) to the conductivity as $\sigma^{s,a;bc}_{(X)}$ for
$X=\mathrm{ii}$, ei, ie, and ee, and we calculate each term. The physical
phenomena are quite different from the case of charge currents if we apply
this reduced density matrix to the calculation of spin currents.

Let us define $\rho_{mn}^{bc\,(X)}{(\omega;\omega_{1},\omega_{2})}$ to satisfy

$\rho_{mn}^{(X)}=\frac{1}{2}\int\frac{d\omega_{1}d\omega_{2}}{(2\pi)^{2}}E^{b}(\omega_{1})E^{c}(\omega_{2})\rho_{mn}^{bc\,(X)}{(\omega;\omega_{1},\omega_{2})}2\pi\delta(\omega-\omega_{1}-\omega_{2})+[(b,\omega_{1})\leftrightarrow(c,\omega_{2})]$
(38)

for $X=2$, ii, ei, ie, and ee. Then,
$\rho_{mn}^{bc\,(2)}{(\omega;\omega_{1},\omega_{2})}$ can be decomposed as
follows Watanabe and Yanase (2021).

$\rho_{mn}^{bc\,(2)}{(\omega;\omega_{1},\omega_{2})}=\sum_{X=\mathrm{ii,\,ei,\,ie,\,ee}}\rho_{mn}^{bc\,(X)}{(\omega;\omega_{1},\omega_{2})}$
(39)
$\rho_{mn}^{bc\,\mathrm{(ii)}}{(\omega;\omega_{1},\omega_{2})}=(-iq)^{2}\delta_{mn}d^{\omega}d^{\omega-\omega_{1}}\partial^{{b}}\partial^{{c}}f_{m}$
(40)
$\rho_{mn}^{bc\,\mathrm{(ei)}}{(\omega;\omega_{1},\omega_{2})}=-iq^{2}d^{\omega}_{mn}d^{\omega-\omega_{1}}\xi^{b}_{mn}\partial^{{c}}f_{mn}$
(41)
$\rho_{mn}^{bc\,\mathrm{(ie)}}{(\omega;\omega_{1},\omega_{2})}=-iq^{2}d^{\omega}_{mn}[\partial^{{b}}(d^{\omega-\omega_{1}}_{mn}f_{mn}\xi^{c}_{mn})-i(\xi^{b}_{mm}-\xi^{b}_{nn})d^{\omega-\omega_{1}}_{mn}f_{mn}\xi^{c}_{mn}]$
(42)
$\rho_{mn}^{bc\,\mathrm{(ee)}}{(\omega;\omega_{1},\omega_{2})}=q^{2}\sum_{\begin{subarray}{c}p\neq
m,n\end{subarray}}d^{\omega}_{mn}[d^{\omega-\omega_{1}}_{pn}\xi^{b}_{mp}\xi^{c}_{pn}f_{np}-d^{\omega-\omega_{1}}_{mp}\xi^{b}_{pn}\xi^{c}_{mp}f_{pm}]$
(43)

Here, we defined

$d^{\omega}_{mn}=\frac{1}{\hbar\omega+i0^{+}-\hbar\omega_{mn}}$ (44)

with $0^{+}$ a positive infinitesimal value that ensures causality of the
response. We write the diagonal components $d^{\omega}_{mm}$ as $d^{\omega}$
since its value does not depend on $m$. Note that the constraint $p\neq m,n$
is included in the definition of $\rho_{mn}^{bc\,\mathrm{(ee)}}$ [Eq. (43)] to
make sure that the intraband effect is excluded. The intraband effect is
considered in $\rho_{mn}^{bc\,\mathrm{(ie)}}$ [Eq. (42)]. Equations (40-43)
corresponds to Eqs. (16-19) of Ref. Watanabe and Yanase (2021).

Using Eq. (34), Eq. (B.1), and Eq. (38) we find the following expression for
the second-order spin photoconductivity.

$\sigma^{s,a;bc}_{(X)}{(\omega;\omega_{1},\omega_{2})}=\frac{q}{2V}\sum_{{\bm{\mathbf{k}}},m,n}j_{mn}^{s,a}\rho_{nm}^{bc\,(X)}{(\omega;\omega_{1},\omega_{2})}+[(b,\omega_{1})\leftrightarrow(c,\omega_{2})].$
(45)

Using the expressions for the second-order density matrix elements [Eqs.
(40-43)], the second-order spin photoconductivity can be computed as follows.

$\sigma^{s,a;bc}{(\omega;\omega_{1},\omega_{2})}=\sigma^{s,a;bc}_{\mathrm{(ii)}}{(\omega;\omega_{1},\omega_{2})}+\sigma^{s,a;bc}_{\mathrm{(ei)}}{(\omega;\omega_{1},\omega_{2})}+\sigma^{s,a;bc}_{\mathrm{(ie)}}{(\omega;\omega_{1},\omega_{2})}+\sigma^{s,a;bc}_{\mathrm{(ee)}}{(\omega;\omega_{1},\omega_{2})}$
(46)
$\sigma^{s,a;bc}_{\mathrm{(ii)}}{(\omega;\omega_{1},\omega_{2})}=\frac{q^{3}}{2V}\sum_{{\bm{\mathbf{k}}},m}-j^{s,a}_{mm}d^{\omega}d^{\omega_{2}}\partial^{{b}}\partial^{{c}}f_{m}+[(b,\omega_{1})\leftrightarrow(c,\omega_{2})]$
(47)
$\sigma^{s,a;bc}_{\mathrm{(ei)}}{(\omega;\omega_{1},\omega_{2})}=\frac{q^{3}}{2V}\sum_{{\bm{\mathbf{k}}},m,n}-ij^{s,a}_{mn}d^{\omega}_{nm}d^{\omega_{2}}\xi^{b}_{nm}\partial^{{c}}f_{nm}+[(b,\omega_{1})\leftrightarrow(c,\omega_{2})]$
(48)
$\sigma^{s,a;bc}_{\mathrm{(ie)}}{(\omega;\omega_{1},\omega_{2})}=\frac{q^{3}}{2V}\sum_{{\bm{\mathbf{k}}},m,n}-ij^{s,a}_{mn}d^{\omega}_{nm}[\partial^{{b}}(d^{\omega_{2}}_{nm}f_{nm}\xi^{c}_{nm})-i(\xi^{b}_{nn}-\xi^{b}_{mm})d^{\omega_{2}}_{nm}f_{nm}\xi^{c}_{nm}]+[(b,\omega_{1})\leftrightarrow(c,\omega_{2})]$
(49)
$\sigma^{s,a;bc}_{\mathrm{(ee)}}{(\omega;\omega_{1},\omega_{2})}=\frac{q^{3}}{2V}\sum_{\begin{subarray}{c}{\bm{\mathbf{k}}},m,n,p\\\
p\neq
m,n\end{subarray}}j^{s,a}_{mn}d^{\omega}_{nm}[d^{\omega_{2}}_{pm}\xi^{b}_{np}\xi^{c}_{pm}f_{mp}-d^{\omega_{2}}_{np}\xi^{b}_{pm}\xi^{c}_{np}f_{pn}]+[(b,\omega_{1})\leftrightarrow(c,\omega_{2})]$
(50)

Now, we take the DC limit [Eq. (36)] and analyze each term.

###  B.2 Drude current

First, we show that the DC limit of
$\sigma^{s,a;bc}_{\mathrm{(ii)}}{(\omega;\omega_{1},\omega_{2})}$ is the Drude
conductivity [Eq. (8)].

$\displaystyle\sigma^{s,a;bc}_{\mathrm{(ii)}}(0;-\Omega,\Omega)=$
$\displaystyle\lim_{\omega\rightarrow
0}\frac{q^{3}}{2\hbar^{2}V\omega}\sum_{{\bm{\mathbf{k}}},m}-j^{s,a}_{mm}\left(\frac{1}{\Omega+\omega}+\frac{1}{-\Omega}\right)\partial^{{b}}\partial^{{c}}f_{m}$
$\displaystyle=$ $\displaystyle\lim_{\omega\rightarrow
0}\frac{q^{3}}{2\hbar^{2}V\omega}\sum_{{\bm{\mathbf{k}}},m}j^{s,a}_{mm}\frac{\omega}{(\Omega+\omega)\Omega}\partial^{{b}}\partial^{{c}}f_{m}$
$\displaystyle=$
$\displaystyle\frac{q^{3}}{2\hbar^{2}V\Omega^{2}}\sum_{{\bm{\mathbf{k}}},m}j^{s,a}_{mm}\partial^{{b}}\partial^{{c}}f_{m}$
$\displaystyle=$ $\displaystyle\sigma^{s,a;bc}_{\mathrm{Drude}}(\Omega)$ (51)

###  B.3 Berry curvature dipole current

Next, we show that the DC limit of
$\sigma^{s,a;bc}_{\mathrm{(ei)}}{(\omega;\omega_{1},\omega_{2})}$ is the Berry
curvature dipole current contribution [Eq. (9)].

$\displaystyle\sigma^{s,a;bc}_{\mathrm{(ei)}}(0;-\Omega,\Omega)=$
$\displaystyle\frac{q^{3}}{2\hbar^{2}V\Omega}\sum_{\begin{subarray}{c}{\bm{\mathbf{k}}},m,n\\\
m\neq
n\end{subarray}}\frac{j^{s,a}_{mn}v^{b}_{nm}}{\omega_{mn}^{2}}\partial^{{c}}f_{nm}+[(b,-\Omega)\leftrightarrow(c,\Omega)]$
$\displaystyle=$
$\displaystyle\frac{q^{3}}{2\hbar^{2}V\Omega}\sum_{\begin{subarray}{c}{\bm{\mathbf{k}}},m,n\\\
m\neq
n\end{subarray}}\left(\frac{j^{s,a}_{nm}v^{b}_{mn}}{\omega_{mn}^{2}}\partial^{{c}}f_{m}-\frac{j^{s,a}_{mn}v^{b}_{nm}}{\omega_{mn}^{2}}\partial^{{c}}f_{m}\right)-(b\leftrightarrow
c)$ $\displaystyle=$
$\displaystyle\frac{-iq^{3}}{\hbar^{2}V\Omega}\sum_{\begin{subarray}{c}{\bm{\mathbf{k}}},m,n\\\
m\neq
n\end{subarray}}\frac{\imaginary(j^{s,a}_{mn}v^{b}_{nm})}{\omega_{mn}^{2}}\partial^{{c}}f_{m}-(b\leftrightarrow
c)$ $\displaystyle=$ $\displaystyle\sigma^{s,a;bc}_{\mathrm{BCD}}(\Omega)$
(52)

Defining the spin Berry curvature Guo  _et al._ (2008)

$\Omega^{s,ab}_{m}=-2\imaginary\sum_{n\neq
m}\frac{j^{s,a}_{mn}v^{b}_{nm}}{\omega_{mn}^{2}},$ (53)

one can rewrite Eq. (B.3) as

$\displaystyle\sigma^{s,a;bc}_{\mathrm{BCD}}=$
$\displaystyle\frac{iq^{3}}{2\hbar^{2}V\Omega}\sum_{{\bm{\mathbf{k}}},m}\Omega^{s,ab}_{m}\partial^{{c}}f_{m}-(b\leftrightarrow
c)=\frac{-iq^{3}}{2\hbar^{2}V\Omega}\sum_{{\bm{\mathbf{k}}},m}f_{m}(\partial^{{c}}\Omega^{s,ab}_{m}-\partial^{{b}}\Omega^{s,ac}_{m}).$
(54)

Equation (54) clearly reveals that $\sigma^{s,a;bc}_{\mathrm{BCD}}$ is
proportional to the momentum space dipole of the spin Berry curvature.

###  B.4 Injection, shift, rectification currents

Finally, we derive the formulas for the injection, shift, and rectification
photoconductivities [Eqs. (10-II.1)] by showing that the sum of
$\sigma^{s,a;bc}_{\mathrm{(ie)}}$ and $\sigma^{s,a;bc}_{\mathrm{(ee)}}$ can be
reorganized as the sum of injection, shift, and rectification currents.

First, let us consider $\sigma^{s,a;bc}_{\mathrm{(ee,d)}}$, which is the
component of $\sigma^{s,a;bc}_{\mathrm{(ee)}}$ [Eq. (50)] with $m=n$:

$\displaystyle\sigma^{s,a;bc}_{\mathrm{(ee,d)}}{(\omega;\omega_{1},\omega_{2})}=$
$\displaystyle\frac{q^{3}}{2V}\sum_{{\bm{\mathbf{k}}},m,p}j^{s,a}_{mm}d^{\omega}[d^{\omega_{2}}_{pm}\xi^{b}_{mp}\xi^{c}_{pm}f_{mp}-d^{\omega_{2}}_{mp}\xi^{b}_{pm}\xi^{c}_{mp}f_{pm}]+[(b,\omega_{1})\leftrightarrow(c,\omega_{2})]$
$\displaystyle=$ $\displaystyle\frac{q^{3}}{2\hbar
V\omega}\sum_{{\bm{\mathbf{k}}},m,p}\Delta
j^{s,a}_{mp}d^{\omega_{2}}_{pm}\xi^{b}_{mp}\xi^{c}_{pm}f_{mp}+[(b,\omega_{1})\leftrightarrow(c,\omega_{2})]$
$\displaystyle=$ $\displaystyle\frac{q^{3}}{2\hbar
V\omega}\sum_{{\bm{\mathbf{k}}},m,p}\Delta
j^{s,a}_{mp}\xi^{b}_{mp}\xi^{c}_{pm}f_{mp}(d^{\omega_{1}}_{mp}+d^{\omega_{2}}_{pm}).$
(55)

Since we are interested in $\sigma^{s,a;bc}(0;-\Omega,\Omega)$, using

$\displaystyle d^{\omega_{1}}_{mp}+d^{\omega_{2}}_{pm}=$
$\displaystyle\frac{1}{-\hbar\Omega+i0^{+}-\hbar\omega_{mp}}+\frac{1}{\hbar(\Omega+\omega)+i0^{+}-\hbar\omega_{pm}}$
$\displaystyle=$ $\displaystyle-2\pi
i\delta(\hbar\Omega+\hbar\omega_{mp})-\frac{\hbar\omega}{(\hbar\Omega+\hbar\omega_{mp})^{2}}+\mathcal{O}(\omega^{2}),$
(56)

we find

$\displaystyle\sigma^{s,a;bc}_{\mathrm{(ee,d)}}(\omega)=$
$\displaystyle\frac{-i\pi
q^{3}}{\hbar^{2}V\omega}\sum_{{\bm{\mathbf{k}}},m,n}\Delta
j^{s,a}_{mn}\xi^{b}_{mn}\xi^{c}_{nm}f_{mn}\delta(\Omega+\omega_{mn})-\frac{q^{3}}{2V}\sum_{{\bm{\mathbf{k}}},m,n}\Delta
j^{s,a}_{mn}\xi^{b}_{mn}\xi^{c}_{nm}f_{mn}\frac{1}{(\hbar\Omega+\hbar\omega_{mn})^{2}}+\mathcal{O}(\omega)$
$\displaystyle=$ $\displaystyle\frac{-i\pi
q^{3}}{\hbar^{2}V\omega\Omega^{2}}\sum_{{\bm{\mathbf{k}}},m,n}f_{mn}\Delta
j^{s,a}_{mn}v^{b}_{mn}v^{c}_{nm}\delta(\Omega+\omega_{mn})-\frac{q^{3}}{2\hbar^{2}V}\sum_{\begin{subarray}{c}{\bm{\mathbf{k}}},m,n\\\
m\neq n\end{subarray}}\frac{f_{mn}}{\omega_{mn}^{2}}\frac{\Delta
j^{s,a}_{mn}v^{b}_{mn}v^{c}_{nm}}{(\Omega+\omega_{mn})^{2}}+\mathcal{O}(\omega).$
(57)

The first term of Eq. (B.4), which contains the $1/\omega$ divergence, is the
injection-current contribution:

$\sigma^{s,a;bc}_{\mathrm{inj.}}(\omega)=\frac{-i\pi
q^{3}}{\hbar^{2}V\omega\Omega^{2}}\sum_{{\bm{\mathbf{k}}},m,n}f_{mn}\Delta
j^{s,a}_{mn}v^{b}_{mn}v^{c}_{nm}\delta(\Omega+\omega_{mn}).$ (58)

The $1/\omega$ divergence in the injection current [Eq. (58)] indicates that
the photocurrent increases linearly with time. This divergence can be avoided
by introducing a phenomenological relaxation time $\tau$ Passos  _et al._
(2018), which is justified by calculations based on the Floquet formalism de
Juan  _et al._ (2017). The DC injection conductivity then becomes

$\displaystyle\sigma^{s,a;bc}_{\mathrm{inj.,\tau}}=\tau\eta^{s,a;bc}_{\mathrm{inj.}}=-\tau\frac{\pi
q^{3}}{\hbar^{2}\Omega^{2}V}\sum_{{\bm{\mathbf{k}}},m,n}f_{mn}\Delta
j^{s,a}_{mn}v^{b}_{mn}v^{c}_{nm}\delta(\Omega+\omega_{mn}).$ (59)

The second term of Eq. (B.4) is a part of the rectification current:

$\sigma^{s,a;bc}_{\mathrm{rect.,1}}=-\frac{q^{3}}{2\hbar^{2}V}\sum_{\begin{subarray}{c}{\bm{\mathbf{k}}},m,n\\\
m\neq n\end{subarray}}\frac{f_{mn}}{\omega_{mn}^{2}}\frac{\Delta
j^{s,a}_{mn}v^{b}_{mn}v^{c}_{nm}}{(\Omega+\omega_{mn})^{2}}$ (60)

Next, let us consider $\sigma^{s,a;bc}_{\mathrm{(ie)}}$. In the DC limit, one
finds

$\displaystyle\sigma^{s,a;bc}_{\mathrm{(ie)}}=$
$\displaystyle\frac{q^{3}}{2V}\sum_{\begin{subarray}{c}{\bm{\mathbf{k}}},m,n\\\
m\neq
n\end{subarray}}-ij^{s,a}_{mn}\frac{1}{\varepsilon_{mn}}[\partial^{b}(d^{\Omega}_{nm}f_{nm}\xi^{c}_{nm})-i(\xi^{b}_{nn}-\xi^{b}_{mm})d^{\Omega}_{nm}f_{nm}\xi^{c}_{nm}]+[(b,-\Omega)\leftrightarrow(c,\Omega)]$
$\displaystyle=$
$\displaystyle\frac{q^{3}}{2V}\sum_{\begin{subarray}{c}{\bm{\mathbf{k}}},m,n\\\
m\neq
n\end{subarray}}i\left[(\partial^{b}j^{s,a}_{mn})-j^{s,a}_{mn}\frac{\Delta^{b}_{mn}}{\omega_{mn}}+ij^{s,a}_{mn}(\xi^{b}_{nn}-\xi^{b}_{mm})\right]\frac{1}{\hbar\omega_{mn}}d^{\Omega}_{nm}f_{nm}\xi^{c}_{nm}+[(b,-\Omega)\leftrightarrow(c,\Omega)]$
$\displaystyle=$ $\displaystyle\frac{iq^{3}}{2\hbar
V}\sum_{\begin{subarray}{c}{\bm{\mathbf{k}}},m,n\\\ m\neq
n\end{subarray}}\left(j^{s,a;b}_{mn}-\frac{j^{s,a}_{mn}\Delta^{b}_{mn}}{\omega_{mn}}\right)\frac{f_{nm}}{\omega_{mn}}d^{\Omega}_{nm}\xi^{c}_{nm}+[(b,-\Omega)\leftrightarrow(c,\Omega)].$
(61)

In the second equality, we used partial integration. We also defined

$\Delta^{b}_{mn}\equiv v^{b}_{mm}-v^{b}_{nn}$ (62)

and the generalized derivative

$j^{s,a;b}_{mn}\equiv\partial^{b}j^{s,a}_{mn}-ij^{s,a}_{mn}(\xi^{b}_{mm}-\xi^{b}_{nn}).$
(63)

Note that the generalized derivative $j^{s,a;b}_{mn}$ is an actual derivative
with respect to the crystal momentum, while the spin-velocity derivative
$d^{s,b;a}_{mn}$ is not.

Last, let us consider the remaining component of
$\sigma^{s,a;bc}_{\mathrm{(ee)}}$ with $m\neq n$, which we denote as
$\sigma^{s,a;bc}_{\mathrm{(ee,od)}}$. In the DC limit,
$\sigma^{s,a;bc}_{\mathrm{(ee,od)}}$ becomes

$\displaystyle\sigma^{s,a;bc}_{\mathrm{(ee,od)}}=$
$\displaystyle\frac{q^{3}}{2V}\sum_{\begin{subarray}{c}{\bm{\mathbf{k}}},m,n,p\\\
m\neq n,p\neq
m,n\end{subarray}}j^{s,a}_{mn}\frac{1}{\hbar\omega_{mn}}(d^{\Omega}_{pm}\xi^{b}_{np}\xi^{c}_{pm}f_{mp}-d^{\Omega}_{np}\xi^{b}_{pm}\xi^{c}_{np}f_{pn})+[(b,-\Omega)\leftrightarrow(c,\Omega)]$
$\displaystyle=$ $\displaystyle\frac{q^{3}}{2\hbar
V}\sum_{\begin{subarray}{c}{\bm{\mathbf{k}}},m,n\end{subarray}}\sum_{\begin{subarray}{c}p\\\
p\neq
m,n\end{subarray}}\left(\frac{j^{s,a}_{mp}\xi^{b}_{pn}}{\omega_{mp}}+\frac{j^{s,a}_{pn}\xi^{b}_{mp}}{\omega_{np}}\right)f_{mn}d^{\Omega}_{nm}\xi^{c}_{nm}+[(b,-\Omega)\leftrightarrow(c,\Omega)].$
(64)

In the second equality, we changed the dummy indices $(m,n,p)$ to $(m,p,n)$
and $(p,n,m)$ for the first and second term in the parentheses, respectively.

To add Eq. (B.4) and Eq. (B.4), we use the following formula:

$\displaystyle\frac{j^{s,a;b}_{mn}}{\omega_{mn}}-j^{s,a}_{mn}\frac{\Delta^{b}_{mn}}{\omega_{mn}^{2}}+i\sum_{\begin{subarray}{c}p\\\
p\neq
m,n\end{subarray}}\left(\frac{j^{s,a}_{mp}\xi^{b}_{pn}}{\omega_{mp}}+\frac{j^{s,a}_{pn}\xi^{b}_{mp}}{\omega_{np}}\right)$
$\displaystyle=$
$\displaystyle\frac{j^{s,ab}_{mn}}{\omega_{mn}}-\frac{i}{\omega_{mn}}\sum_{p\neq
n}j^{s,a}_{mp}\xi^{b}_{pn}+\frac{i}{\omega_{mn}}\sum_{p\neq
m}\xi^{b}_{mp}j^{s,a}_{pn}-j^{s,a}_{mn}\frac{\Delta^{b}_{mn}}{\omega_{mn}^{2}}+i\sum_{p\neq
m,n}\frac{j^{s,a}_{mp}\xi^{b}_{pn}}{\omega_{mp}}+i\sum_{p\neq
m,n}\frac{j^{s,a}_{pn}\xi^{b}_{mp}}{\omega_{np}}$ $\displaystyle=$
$\displaystyle\frac{j^{s,ab}_{mn}}{\omega_{mn}}-i\frac{\Delta
j^{s,a}_{mn}\xi^{b}_{mn}}{\omega_{mn}}-\frac{j^{s,a}_{mn}\Delta^{b}_{mn}}{\omega_{mn}^{2}}+\frac{i}{\omega_{mn}}\sum_{p\neq
m,n}\frac{j^{s,a}_{mp}\xi^{b}_{pn}\omega_{pn}}{\omega_{mp}}+\frac{i}{\omega_{mn}}\sum_{p\neq
m,n}\frac{j^{s,a}_{pn}\xi^{b}_{mp}\omega_{mp}}{\omega_{np}}$ $\displaystyle=$
$\displaystyle\frac{j^{s,ab}_{mn}}{\omega_{mn}}-\frac{\Delta
j^{s,a}_{mn}v^{b}_{mn}}{\omega_{mn}^{2}}-\frac{j^{s,a}_{mn}\Delta^{b}_{mn}}{\omega_{mn}^{2}}+\frac{1}{\omega_{mn}}\sum_{p\neq
m,n}\frac{j^{s,a}_{mp}v^{b}_{pn}}{\omega_{mp}}+\frac{1}{\omega_{mn}}\sum_{p\neq
m,n}\frac{v^{b}_{mp}j^{s,a}_{pn}}{\omega_{np}}$ $\displaystyle=$
$\displaystyle\frac{1}{\omega_{mn}}\left(j^{s,ab}_{mn}+\sum_{p\neq
m}\frac{j^{s,a}_{mp}v^{b}_{pn}}{\omega_{mp}}+\sum_{p\neq
n}\frac{v^{b}_{mp}j^{s,a}_{pn}}{\omega_{np}}\right)-\frac{\Delta
j^{s,a}_{mn}v^{b}_{mn}}{\omega_{mn}^{2}}$ $\displaystyle=$
$\displaystyle\frac{d^{s,b;a}_{mn}}{\omega_{mn}}-\frac{\Delta
j^{s,a}_{mn}v^{b}_{mn}}{\omega_{mn}^{2}}.$ (65)

From Eqs. (B.4, B.4, B.4), we find

$\displaystyle\sigma^{s,a;bc}_{\mathrm{(ie)}}+\sigma^{s,a;bc}_{\mathrm{(ee,od)}}$
$\displaystyle=$ $\displaystyle\frac{iq^{3}}{2\hbar
V}\sum_{\begin{subarray}{c}{\bm{\mathbf{k}}},m,n\\\ m\neq
n\end{subarray}}\left[\frac{j^{s,a;b}_{mn}}{\omega_{mn}}-j^{s,a}_{mn}\frac{\Delta^{b}_{mn}}{\omega_{mn}^{2}}+i\sum_{\begin{subarray}{c}c\\\
c\neq
a,b\end{subarray}}\left(\frac{j^{s,a}_{mp}\xi^{b}_{pn}}{\omega_{mp}}+\frac{j^{s,a}_{pn}\xi^{b}_{mp}}{\omega_{np}}\right)\right]f_{nm}d^{\Omega}_{nm}\xi^{c}_{nm}+[(b,-\Omega)\leftrightarrow(c,\Omega)]$
$\displaystyle=$ $\displaystyle\frac{iq^{3}}{2\hbar
V}\sum_{\begin{subarray}{c}{\bm{\mathbf{k}}},m,n\\\ m\neq
n\end{subarray}}\left(\frac{d^{s,b;a}_{mn}}{\omega_{mn}}-\frac{\Delta
j^{s,a}_{mn}v^{b}_{mn}}{\omega_{mn}^{2}}\right)f_{nm}d^{\Omega}_{nm}\xi^{c}_{nm}+[(b,-\Omega)\leftrightarrow(c,\Omega)].$
(66)

Now, let us separate $d^{\Omega}_{nm}$ into the delta function part and the
principal value part:

$d^{\Omega}_{nm}=\frac{1}{\hbar\Omega+i0^{+}+\hbar\omega_{mn}}=-i\pi\delta(\hbar\Omega+\hbar\omega_{mn})+\mathrm{P}\frac{1}{\hbar(\Omega+\omega_{mn})}.$
(67)

The delta function part is the shift current:

$\displaystyle\sigma^{s,a;bc}_{\mathrm{shift}}=$ $\displaystyle\frac{i\pi
q^{3}}{2\hbar^{2}V}\sum_{\begin{subarray}{c}{\bm{\mathbf{k}}},m,n\\\ m\neq
n\end{subarray}}\left(\frac{d^{s,b;a}_{mn}v^{c}_{nm}}{\omega_{mn}^{2}}-\frac{\Delta
j^{s,a}_{mn}v^{b}_{mn}v^{c}_{nm}}{\omega_{mn}^{3}}\right)f_{nm}\delta(\Omega+\omega_{mn})+[(b,-\Omega)\leftrightarrow(c,\Omega)]$
$\displaystyle=$ $\displaystyle\frac{i\pi
q^{3}}{2\hbar^{2}V\Omega^{2}}\sum_{{\bm{\mathbf{k}}},m,n}(d^{s,b;a}_{mn}v^{c}_{nm}-d^{s,c;a}_{nm}v^{b}_{mn})f_{nm}\delta(\Omega+\omega_{mn}).$
(68)

The principal value part is

$\displaystyle\sigma^{s,a;bc}_{\mathrm{rect.,2}}=$
$\displaystyle\frac{q^{3}}{2\hbar^{2}V}\sum_{\begin{subarray}{c}{\bm{\mathbf{k}}},m,n\\\
m\neq
n\end{subarray}}f_{mn}\left(\frac{d^{s,b;a}_{mn}v^{c}_{nm}}{\omega_{mn}^{2}}-\frac{\Delta
j^{s,a}_{mn}v^{b}_{mn}v^{c}_{nm}}{\omega_{mn}^{3}}\right)\mathrm{P}\frac{1}{\Omega+\omega_{mn}}+[(b,-\Omega)\leftrightarrow(c,\Omega)]$
$\displaystyle=$
$\displaystyle\frac{q^{3}}{2\hbar^{2}V}\sum_{\begin{subarray}{c}{\bm{\mathbf{k}}},m,n\\\
m\neq
n\end{subarray}}f_{mn}\left(\frac{d^{s,b;a}_{mn}v^{c}_{nm}+d^{s,c;a}_{nm}v^{b}_{mn}}{\omega_{mn}^{2}}-\frac{2\Delta
j^{s,a}_{mn}v^{b}_{mn}v^{c}_{nm}}{\omega_{mn}^{3}}\right)\mathrm{P}\frac{1}{\Omega+\omega_{mn}}$
(69)

Adding Eq. (60) and Eq. (B.4), we obtain the rectification current:

$\displaystyle\sigma^{s,a;bc}_{\mathrm{rect.}}=$
$\displaystyle\sigma^{s,a;bc}_{\mathrm{rect.,1}}+\sigma^{s,a;bc}_{\mathrm{rect.,2}}$
$\displaystyle=$
$\displaystyle\frac{q^{3}}{2\hbar^{2}V}\sum_{\begin{subarray}{c}{\bm{\mathbf{k}}},m,n\\\
m\neq
n\end{subarray}}f_{mn}\left[\left(d^{s,b;a}_{mn}v^{c}_{nm}+d^{s,c;a}_{nm}v^{b}_{mn}-\frac{2\Delta
j^{s,a}_{mn}v^{b}_{mn}v^{c}_{nm}}{\omega_{mn}}\right)\frac{1}{\omega_{mn}^{2}}\mathrm{P}\frac{1}{\Omega+\omega_{mn}}-\frac{\Delta
j^{s,a}_{mn}v^{b}_{mn}v^{c}_{nm}}{\omega_{mn}^{2}(\Omega+\omega_{mn})^{2}}\right].$
(70)

##  Appendix C Symmetry analysis of photoconductivity

In this section, we study the symmetry properties of the second-order
photoconductivity tensors. The real symmetric (imaginary anti-symmetric)
component of the photoconductivity tensor,
$\real(\sigma^{s,a;bc}+\sigma^{s,a;cb})$
($\imaginary(\sigma^{s,a;bc}-\sigma^{s,a;cb})$) corresponds to the generation
of a real-valued current under linearly (circularly) polarized light. In other
words, the photoconductivity tensor for linearly (circularly) polarized light
is symmetric (anti-symmetric) with respect to the exchange of indices $b$ and
$c$. For later use, we define

$\sigma^{s,a;bc}_{\rm
L}=\frac{1}{2}\real\left(\sigma^{s,a;bc}+\sigma^{s,a;cb}\right)$ (71)

and

$\sigma^{s,a;bc}_{\rm
C}=\frac{1}{2}\imaginary\left(\sigma^{s,a;bc}-\sigma^{s,a;cb}\right)\,.$ (72)

The symmetry of the charge shift and injection conductivities were studied in
Ref. Ahn  _et al._ (2020). Here, we complete the analysis by studying the
symmetry of the Drude, BCD, and rectification conductivities as well as the
spin photoconductivities. Let us begin with the charge case.

Following Ref. Ahn  _et al._ (2020), we first study how the matrix elements
transform under a point-group symmetry operation and the time-reversal
operation. We do not assume that the system is invariant under the symmetry
operations. Instead, we study the relationship between the matrix element and
conductivity tensors of the transformed and the original systems. We denote
the quantities of the transformed system with a prime.

First, let us consider a point-group symmetry operation ${\hat{\mathcal{M}}}$
whose real $3\times 3$ transformation matrix is ${\mathcal{M}}$. We denote the
Hamiltonian of the original and transformed systems by $H_{0}$ and
$H_{0}^{\prime}={\hat{\mathcal{M}}}\,H_{0}\,{\hat{\mathcal{M}}}^{-1}$,
respectively. Also, $H_{0}({\bf k})=\exp(-i{\bf k}\cdot{\bf
r})\,H_{0}\,\exp(i{\bf k}\cdot{\bf r})$ and $H_{0}^{\prime}({\bf
k})=\exp(-i{\bf k}\cdot{\bf r})\,H_{0}^{\prime}\,\exp(i{\bf k}\cdot{\bf r})$.
The periodic parts of the Bloch states of the two systems satisfy
$\ket{u^{\prime}_{n\mathbf{k}}}={\hat{\mathcal{M}}}\ket{u_{n{\mathcal{M}}^{-1}{\bm{\mathbf{k}}}}}$.
The eigenvalues satisfy
$\varepsilon^{\prime}_{n\mathbf{k}}=\varepsilon_{n{\mathcal{M}}^{-1}{\bm{\mathbf{k}}}}$.
The velocity matrix element transforms as a vector:

$\displaystyle v_{mn}^{\prime a}({\bm{\mathbf{k}}})=$
$\displaystyle\hbar^{-1}\matrixelement{u^{\prime}_{m\mathbf{k}}}{\partial^{{a}}H_{0}^{\prime}({\bm{\mathbf{k}}})}{u^{\prime}_{n\mathbf{k}}}$
$\displaystyle=$
$\displaystyle\hbar^{-1}\matrixelement{{\hat{\mathcal{M}}}u_{m{\mathcal{M}}^{-1}{\bm{\mathbf{k}}}}}{\partial^{{a}}\left[\exp(-i{\bf
k}\cdot{\bf
r})\,{\hat{\mathcal{M}}}\,H_{0}\,{\hat{\mathcal{M}}}^{-1}\,\exp(i{\bf
k}\cdot{\bf
r})\,\right]}{{\hat{\mathcal{M}}}u_{n{\mathcal{M}}^{-1}{\bm{\mathbf{k}}}}}$
$\displaystyle=$
$\displaystyle\hbar^{-1}\matrixelement{{\hat{\mathcal{M}}}u_{m{\mathcal{M}}^{-1}{\bm{\mathbf{k}}}}}{{\hat{\mathcal{M}}}\,\partial^{{a}}\left[\exp(-i{\bf
k}\cdot{\mathcal{M}}{\bf r})\,H_{0}\,\exp(i{\bf k}\cdot{\mathcal{M}}{\bf
r})\,\right]\,{\hat{\mathcal{M}}}^{-1}}{{\hat{\mathcal{M}}}u_{n{\mathcal{M}}^{-1}{\bm{\mathbf{k}}}}}$
$\displaystyle=$
$\displaystyle\hbar^{-1}\matrixelement{{\hat{\mathcal{M}}}u_{m{\mathcal{M}}^{-1}{\bm{\mathbf{k}}}}}{{\hat{\mathcal{M}}}\,\partial^{{a}}\left[H_{0}({\mathcal{M}}^{-1}{\bm{\mathbf{k}}})\right]\,{\hat{\mathcal{M}}}^{-1}}{{\hat{\mathcal{M}}}u_{n{\mathcal{M}}^{-1}{\bm{\mathbf{k}}}}}$
$\displaystyle=$
$\displaystyle\frac{\partial({\mathcal{M}}^{-1}{\bm{\mathbf{k}}})_{b}}{\partial
k_{a}}\,\hbar^{-1}\matrixelement{u_{m{\mathcal{M}}^{-1}{\bm{\mathbf{k}}}}}{(\partial^{{b}}H_{0})({\mathcal{M}}^{-1}{\bm{\mathbf{k}}})}{u_{n{\mathcal{M}}^{-1}{\bm{\mathbf{k}}}}}$
$\displaystyle=$
$\displaystyle{{\mathcal{M}}^{-1}}_{ba}\,v^{b}_{mn}\left({\mathcal{M}}^{-1}{\bm{\mathbf{k}}}\right)$
$\displaystyle=$
$\displaystyle{\mathcal{M}}_{ab}\,v^{b}_{mn}({\mathcal{M}}^{-1}{\bm{\mathbf{k}}})\,.$
(73)

The derivative operation also transforms as a vector:

$(\partial^{{a}}f)^{\prime}({\bm{\mathbf{k}}})=\partial^{{a}}f({\mathcal{M}}^{-1}{\bm{\mathbf{k}}})=\frac{\partial({\mathcal{M}}^{-1}{\bm{\mathbf{k}}})_{b}}{\partial
k_{a}}(\partial^{{b}}f)({\mathcal{M}}^{-1}{\bm{\mathbf{k}}})={\mathcal{M}}_{ab}(\partial^{{b}}f)({\mathcal{M}}^{-1}{\bm{\mathbf{k}}}).$
(74)

Similarly, one can easily show that the generalized derivative
$d^{0,b;a}_{mn}$ transforms as a rank-2 tensor,

$\displaystyle
d^{\,\prime\,0,b^{\prime};a^{\prime}}_{mn}({\bm{\mathbf{k}}})={\mathcal{M}}_{a^{\prime}a}{\mathcal{M}}_{b^{\prime}b}d^{0,b;a}_{mn}({\mathcal{M}}^{-1}{\bm{\mathbf{k}}})\,,$
(75)

and that each of all the charge photoconductivity tensors (Drude, BCD, shift,
injection, and rectification conductivity tensors) transforms under
${\hat{\mathcal{M}}}$ as

$\sigma^{\,\prime\,0,a^{\prime};b^{\prime}c^{\prime}}={\mathcal{M}}_{a^{\prime}a}{\mathcal{M}}_{b^{\prime}b}{\mathcal{M}}_{c^{\prime}c}\sigma^{0,{a;bc}}\,.$
(76)

Thus, the second-order charge photoconductivity $\sigma^{0,a;bc}$ is a rank-3
tensor.

Next, let us consider the time-reversal operation ${\hat{\mathcal{T}}}$. The
periodic parts of the Bloch states of the transformed and the original systems
satisfy

$\ket{u^{\prime}_{n\mathbf{k}}}={\hat{\mathcal{T}}}\ket{u_{n-{\bm{\mathbf{k}}}}}.$
(77)

The velocity matrix element transforms as

$\displaystyle v_{mn}^{\prime a}({\bm{\mathbf{k}}})=$
$\displaystyle\hbar^{-1}\matrixelement{u^{\prime}_{m\mathbf{k}}}{\partial^{{a}}H_{0}^{\prime}({\bm{\mathbf{k}}})}{u^{\prime}_{n\mathbf{k}}}$
$\displaystyle=$
$\displaystyle\hbar^{-1}\matrixelement{{\hat{\mathcal{T}}}u_{m-{\bm{\mathbf{k}}}}}{\partial^{{a}}\left[{\hat{\mathcal{T}}}H_{0}(-{\bm{\mathbf{k}}}){\hat{\mathcal{T}}}^{-1}\right]}{{\hat{\mathcal{T}}}u_{n-{\bm{\mathbf{k}}}}}$
$\displaystyle=$
$\displaystyle\hbar^{-1}\innerproduct{{\hat{\mathcal{T}}}u_{m-{\bm{\mathbf{k}}}}}{{\hat{\mathcal{T}}}u_{m-{\bm{\mathbf{k}}}}}{{\hat{\mathcal{T}}}\partial^{{a}}\left[H_{0}(-{\bm{\mathbf{k}}})\right]u_{n-{\bm{\mathbf{k}}}}}$
$\displaystyle=$
$\displaystyle\hbar^{-1}\innerproduct{\partial^{{a}}\left[H_{0}(-{\bm{\mathbf{k}}})\right]u_{n-{\bm{\mathbf{k}}}}}{\partial^{{a}}\left[H_{0}(-{\bm{\mathbf{k}}})\right]u_{n-{\bm{\mathbf{k}}}}}{u_{m-{\bm{\mathbf{k}}}}}$
$\displaystyle=$
$\displaystyle-\hbar^{-1}\innerproduct{(\partial^{{a}}H_{0})(-{\bm{\mathbf{k}}})u_{n-{\bm{\mathbf{k}}}}}{(\partial^{{a}}H_{0})(-{\bm{\mathbf{k}}})u_{n-{\bm{\mathbf{k}}}}}{u_{m-{\bm{\mathbf{k}}}}}$
$\displaystyle=$ $\displaystyle-v_{nm}^{a}(-{\bm{\mathbf{k}}})$ (78)

In the fourth equality, we used
$\innerproduct{{\hat{\mathcal{T}}}u}{{\hat{\mathcal{T}}}v}=\innerproduct{v}{u}$.
The derivative operation also transforms as

$(\partial^{{a}}f)^{\prime}({\bm{\mathbf{k}}})=\partial^{{a}}f(-{\bm{\mathbf{k}}})=-(\partial^{{a}}f)(-{\bm{\mathbf{k}}})\,.$
(79)

The generalized derivative $d^{0,b;a}_{mn}$ transforms as

$\displaystyle d^{\,\prime\,0,b;a}_{mn}({\bm{\mathbf{k}}})=$ $\displaystyle
j^{0,ab}_{nm}(-{\bm{\mathbf{k}}})$ $\displaystyle+\sum_{p\neq
m}\frac{v^{a}_{pm}(-{\bm{\mathbf{k}}})v^{b}_{np}(-{\bm{\mathbf{k}}})}{\omega_{mp}(-{\bm{\mathbf{k}}})}$
$\displaystyle+\sum_{p\neq
n}\frac{v^{b}_{pm}(-{\bm{\mathbf{k}}})v^{a}_{np}(-{\bm{\mathbf{k}}})}{\omega_{np}(-{\bm{\mathbf{k}}})}$
$\displaystyle=$ $\displaystyle d^{0,b;a}_{nm}(-{\bm{\mathbf{k}}}).$ (80)

Now, using the symmetry properties of the matrix elements, we study the
transformation of charge conductivity tensors under ${\hat{\mathcal{T}}}$.

###  C.1 Drude

Equation (8) shows that the charge Drude conductivity is always symmetric in
the exchange of $b$ and $c$. Hence, we find

$\sigma^{0,{a;bc}}_{\rm
Drude,C}=\frac{1}{2}\imaginary\left(\sigma^{0,{a;bc}}_{\rm
Drude}-\sigma^{0,a;cb}_{\rm Drude}\right)=0.$ (81)

Further from Eq. (8), one can easily show that the charge Drude conductivity
transforms as

$\sigma^{\,\prime\,0,{a;bc}}_{\mathrm{Drude}}=-\sigma^{0,{a;bc}}_{\mathrm{Drude}}.$
(82)

The $(-1)$ factor shows that the charge Drude conductivity is odd under
${\hat{\mathcal{T}}}$. Thus, in a ${\hat{\mathcal{T}}}$-symmetric system, the
charge Drude conductivity is zero because
$\sigma^{0,a;bc}_{\mathrm{Drude}}=\sigma^{\,\prime\,0,a;bc}_{\mathrm{Drude}}=-\sigma^{0,a;bc}_{\mathrm{Drude}}$
holds.

One can use $j^{0,a}_{mm}=v^{a}_{mm}=\partial^{{a}}\varepsilon_{m}/\hbar$ and
partial integration to rewrite the charge Drude conductivity as follows:

$\displaystyle\sigma^{0,{a;bc}}_{\mathrm{Drude}}=\frac{q^{3}}{2\hbar^{2}V\Omega^{2}}\sum_{{\bm{\mathbf{k}}},m}j^{0,a}_{mm}\partial^{{b}}\partial^{{c}}f_{m}=\frac{q^{3}}{2\hbar^{3}V\Omega^{2}}\sum_{{\bm{\mathbf{k}}},m}(\partial^{{a}}\varepsilon_{m})(\partial^{{b}}\partial^{{c}}f_{m})=\frac{q^{3}}{2\hbar^{3}V\Omega^{2}}\sum_{{\bm{\mathbf{k}}},m}(\partial^{{a}}\partial^{{b}}\partial^{{c}}\varepsilon_{m})f_{m}.$
(83)

Hence, the charge Drude current is symmetric under the permutation of all
three indices.

###  C.2 Berry curvature dipole

For the BCD current, Eq. (9) shows that it is always anti-symmetric in the
exchange of $b$ and $c$. Hence, we find

$\sigma^{0,{a;bc}}_{\rm BCD,L}=\frac{1}{2}\real\left(\sigma^{0,{a;bc}}_{\rm
BCD}+\sigma^{0,a;cb}_{\rm BCD}\right)=0.$ (84)

The summand of Eq. (9) transforms under ${\hat{\mathcal{T}}}$ as

$\left[\frac{\imaginary(j_{mn}^{0,a}v_{nm}^{b})}{\omega_{mn}^{2}}(\partial^{{c}}f_{m})\right]^{\prime}({\bm{\mathbf{k}}})=\left[\frac{\imaginary(j_{nm}^{0,a}v_{mn}^{b})}{\omega_{mn}^{2}}(-\partial^{{c}}f_{m})\right](-{\bm{\mathbf{k}}})=\left[\frac{\imaginary(j_{mn}^{0,a}v_{nm}^{b})}{\omega_{mn}^{2}}(\partial^{{c}}f_{m})\right](-{\bm{\mathbf{k}}})$
(85)

Hereafter, we use the notation that the arguments ${\bm{\mathbf{k}}}$ and
$-{\bm{\mathbf{k}}}$ as well as the prime apply to all matrix elements inside
the outermost parentheses on its left. In the last equality of Eq. (85), we
used the Hermiticity of the current operators. Thus, the charge BCD
conductivity is symmetric under ${\hat{\mathcal{T}}}$:

$\sigma^{\,\prime\,0,{a;bc}}_{\mathrm{BCD}}=\sigma^{0,{a;bc}}_{\mathrm{BCD}}.$
(86)

###  C.3 Shift

Next, let us consider the shift conductivity. Using Eqs. (C, C), we find that
the matrix element in Eq. (II.1) transforms as

$\left(d^{0,b;a}_{mn}v^{c}_{nm}-d^{0,c;a}_{nm}v^{b}_{mn}\right)^{\prime}({\bm{\mathbf{k}}})=-\left(d^{0,b;a}_{nm}v^{c}_{mn}-d^{0,c;a}_{mn}v^{b}_{nm}\right)(-{\bm{\mathbf{k}}})=\left(d^{0,c;a}_{mn}v^{b}_{nm}-d^{0,b;a}_{nm}v^{c}_{mn}\right)(-{\bm{\mathbf{k}}}).$
(87)

Note that the time reversal exchanges the field direction indices $b$ and $c$
without a sign flip. Thus, the charge shift conductivity transforms under
${\hat{\mathcal{T}}}$ as

$\sigma^{\,\prime\,0,{a;cb}}_{\mathrm{shift}}=\sigma^{0,{a;bc}}_{\mathrm{shift}}.$
(88)

The shift conductivity for linearly polarized light transforms under
${\hat{\mathcal{T}}}$ as

$\sigma^{\,\prime\,0,{a;bc}}_{\mathrm{shift,L}}=\frac{1}{2}\real\left(\sigma^{\,\prime\,0,{a;bc}}_{\mathrm{shift}}+\sigma^{\,\prime\,0,{a;cb}}_{\mathrm{shift}}\right)=\frac{1}{2}\real\left(\sigma^{0,{a;cb}}_{\mathrm{shift}}+\sigma^{0,{a;bc}}_{\mathrm{shift}}\right)=\sigma^{0,{a;bc}}_{\mathrm{shift,L}},$
(89)

while that for circularly polarized light transforms as

$\sigma^{\,\prime\,0,{a;bc}}_{\mathrm{shift,C}}=\frac{1}{2}\imaginary\left(\sigma^{\,\prime\,0,{a;bc}}_{\mathrm{shift}}-\sigma^{\,\prime\,0,{a;cb}}_{\mathrm{shift}}\right)=\frac{1}{2}\imaginary\left(\sigma^{0,{a;cb}}_{\mathrm{shift}}-\sigma^{0,{a;bc}}_{\mathrm{shift}}\right)=-\sigma^{0,{a;bc}}_{\mathrm{shift,C}}.$
(90)

The shift conductivity for linearly (circularly) polarized light is even (odd)
under time reversal.

###  C.4 Injection

Next, let us consider the injection conductivity. Using Eqs. (C), we find that
the matrix elements in Eq. (10) transforms under ${\hat{\mathcal{T}}}$ as

$\left(f_{mn}\Delta
j^{0,a}_{mn}v^{b}_{mn}v^{c}_{nm}\right)^{\prime}({\bm{\mathbf{k}}})=-\left(f_{mn}\Delta
j^{0,a}_{mn}v^{b}_{nm}v^{c}_{mn}\right)(-{\bm{\mathbf{k}}})$ (91)

The time reversal exchanges the field direction indices $b$ and $c$ and gives
a sign flip. Thus, the rectification conductivity transforms under
${\hat{\mathcal{T}}}$ as

$\sigma^{\,\prime\,0,{a;cb}}_{\mathrm{inj.}}=-\sigma^{0,{a;bc}}_{\mathrm{inj.}}.$
(92)

The injection conductivity for linearly polarized light transforms under
${\hat{\mathcal{T}}}$ as

$\sigma^{\,\prime\,0,{a;bc}}_{\mathrm{inj.,L}}=\frac{1}{2}\real\left(\sigma^{\,\prime\,0,{a;bc}}_{\mathrm{inj.}}+\sigma^{\,\prime\,0,{a;cb}}_{\mathrm{inj.}}\right)=-\frac{1}{2}\real\left(\sigma^{0,{a;cb}}_{\mathrm{inj.}}+\sigma^{0,{a;bc}}_{\mathrm{inj.}}\right)=-\sigma^{0,{a;bc}}_{\mathrm{inj.,L}},$
(93)

while that for circularly polarized light transforms as

$\sigma^{\,\prime\,0,{a;bc}}_{\mathrm{inj.,C}}=\frac{1}{2}\imaginary\left(\sigma^{\,\prime\,0,{a;bc}}_{\mathrm{inj.}}-\sigma^{\,\prime\,0,{a;cb}}_{\mathrm{inj.}}\right)=-\frac{1}{2}\imaginary\left(\sigma^{0,{a;cb}}_{\mathrm{inj.}}-\sigma^{0,{a;bc}}_{\mathrm{inj.}}\right)=\sigma^{0,{a;bc}}_{\mathrm{inj.,C}}.$
(94)

Note the sign difference between Eqs. (89, 90) and Eqs. (93, 94). The
injection conductivity for linearly (circularly) polarized light is odd (even)
under time reversal.

###  C.5 Rectification

Finally, let us consider the rectification conductivity. Using Eqs. (C, C),
one can show the term in the inner parentheses of Eq. (II.1) transforms under
${\hat{\mathcal{T}}}$ as

$\left(d^{0,b;a}_{mn}v^{c}_{nm}+d^{0,c;a}_{nm}v^{b}_{mn}-\frac{2\Delta
j^{0,a}_{mn}v^{b}_{mn}v^{c}_{nm}}{\omega_{mn}}\right)^{\prime}({\bm{\mathbf{k}}})=-\left(d^{0,b;a}_{nm}v^{c}_{mn}+d^{0,c;a}_{mn}v^{b}_{nm}-\frac{2\Delta
j^{0,a}_{mn}v^{b}_{nm}v^{c}_{mn}}{\omega_{mn}}\right)(-{\bm{\mathbf{k}}}).$
(95)

The time reversal exchanges the field direction indices $b$ and $c$ and gives
a sign flip. Thus, the rectification conductivity transforms under
${\hat{\mathcal{T}}}$ as

$\sigma^{\,\prime\,0,{a;cb}}_{\mathrm{rect.}}=-\sigma^{0,{a;bc}}_{\mathrm{rect.}},$
(96)

which has the same form as Eq. (92). Thus, the symmetry property of the
rectification conductivity is identical to that of the injection conductivity.
This result completes the symmetry analysis of second-order charge
photocurrents.

###  C.6 Spin shift vector

Now, let us discuss the symmetry transformation property of the spin shift
vector, which we have proposed as a key quantity that describes the spin shift
and rectification currents. Note that in the shift vector
${\widetilde{R}}^{s,{\mathcal{E}};a}_{mn}$, one should not treat the light
polarization ${\bm{\mathbf{{\mathcal{E}}}}}$ as a tensor index because the
shift vector does not transform like a tensor for that index (see the second
term of Eq. (20)). Hence, we treat ${\bm{\mathbf{{\mathcal{E}}}}}$ as an
external parameter that is transformed together with the system.

First, let us consider a rotation operation $\hat{{\mathcal{M}}}$, whose real
$3\times 3$ transformation matrix is ${\mathcal{M}}$. The light polarization
vector transforms as
${\mathcal{E}}^{{}^{\prime}a}={\mathcal{M}}_{ab}{\mathcal{E}}^{b}$. If an
operator $\mathcal{O}^{a}({\bm{\mathbf{k}}})$ transforms as
$\mathcal{O}^{\prime
a}_{mn}({\bm{\mathbf{k}}})={\mathcal{M}}_{ab}\mathcal{O}^{b}_{mn}({\mathcal{M}}^{-1}{\bm{\mathbf{k}}})$
under a rotation ${\mathcal{M}}$, we find that
$\mathcal{O}^{{\mathcal{E}}}_{mn}$ transforms as

$\left(\mathcal{O}^{{\mathcal{E}}}_{mn}\right)^{\prime}({\bm{\mathbf{k}}})=\mathcal{O}^{\prime
a}_{mn}({\bm{\mathbf{k}}}){\mathcal{E}}^{\prime
a}={\mathcal{M}}_{ab}{\mathcal{M}}_{ac}\mathcal{O}^{b}_{mn}({\mathcal{M}}^{-1}{\bm{\mathbf{k}}}){\mathcal{E}}^{c}=\delta_{b,c}\mathcal{O}^{b}_{mn}({\mathcal{M}}^{-1}{\bm{\mathbf{k}}}){\mathcal{E}}^{c}=\mathcal{O}^{{\mathcal{E}}}_{mn}({\mathcal{M}}^{-1}{\bm{\mathbf{k}}}).$
(97)

In other words, the inner product of two vector quantities
[$\mathcal{O}^{a}_{mn}({\bm{\mathbf{k}}})$ and ${\mathcal{E}}^{a}$] transforms
as a scalar. Then, one finds that the spin shift vector [Eq. (20)] transforms
as

$\left({\widetilde{R}}^{s^{\prime},{\mathcal{E}};a^{\prime}}_{mn}\right)^{\prime}({\bm{\mathbf{k}}})=i\left(\frac{d^{s^{\prime},{\mathcal{E}};a^{\prime}}_{mn}}{v^{{\mathcal{E}}}_{mn}}-\frac{\Delta
j^{s^{\prime},a^{\prime}}_{mn}}{\omega_{mn}}\right)^{\prime}({\bm{\mathbf{k}}})=i{\mathcal{M}}_{s^{\prime}s}{\mathcal{M}}_{a^{\prime}a}\left(\frac{d^{s,{\mathcal{E}};a}_{mn}}{v^{\mathcal{E}}_{mn}}-\frac{\Delta
j^{s,a}_{mn}}{\omega_{mn}}\right)({\mathcal{M}}^{-1}{\bm{\mathbf{k}}})={\mathcal{M}}_{s^{\prime}s}{\mathcal{M}}_{a^{\prime}a}{\widetilde{R}}^{s,{\mathcal{E}};a}_{mn}({\mathcal{M}}^{-1}{\bm{\mathbf{k}}}).$
(98)

For an improper rotation, the spin shift vector gets an additional $-1$ factor
because the spin-current operator is an axial vector.

Next, for the time-reversal operation, the light polarization transforms as
${\mathcal{E}}^{\prime a}=({\mathcal{E}}^{a})^{*}$. If $\mathcal{O}^{a}$
transforms under time reversal as $\mathcal{O}^{\prime
a}_{mn}({\bm{\mathbf{k}}})=\pm\left[\mathcal{O}^{a}_{mn}(-{\bm{\mathbf{k}}})\right]^{*}$,
$\mathcal{O}^{\mathcal{E}}$ transforms as

$\left(\mathcal{O}^{{\mathcal{E}}}_{mn}\right)^{\prime}({\bm{\mathbf{k}}})=\mathcal{O}^{\prime
a}_{mn}({\bm{\mathbf{k}}}){\mathcal{E}}^{\prime
a}=\pm\left[\mathcal{O}^{a}_{mn}(-{\bm{\mathbf{k}}})\right]^{*}({\mathcal{E}}^{a})^{*}=\pm\left[\mathcal{O}^{{\mathcal{E}}}_{mn}(-{\bm{\mathbf{k}}})\right]^{*}.$
(99)

Then, for $s\neq 0$ we find

$\left({\widetilde{R}}^{s,{\mathcal{E}};a}_{mn}\right)^{\prime}({\bm{\mathbf{k}}})=\left[i\left(\frac{d^{s,{\mathcal{E}};a}_{mn}}{v^{{\mathcal{E}}}_{mn}}-\frac{\Delta
j^{s,a}_{mn}}{\omega_{mn}}\right)\right]^{\prime}({\bm{\mathbf{k}}})=i\left[\left(\frac{d^{s,{\mathcal{E}};a}_{mn}}{v^{{\mathcal{E}}}_{mn}}-\frac{\Delta
j^{s,a}_{mn}}{\omega_{mn}}\right)(-{\bm{\mathbf{k}}})\right]^{*}=-\left({\widetilde{R}}^{s,{\mathcal{E}};a}_{mn}(-{\bm{\mathbf{k}}})\right)^{*}.$
(100)

In the second equality, we included the $-1$ factor that occurs because the
spin-current operator is odd under time reversal. For the charge case, one
finds

$\left({\widetilde{R}}^{0,{\mathcal{E}}^{;}a}_{mn}\right)^{\prime}({\bm{\mathbf{k}}})=\left({\widetilde{R}}^{0,{\mathcal{E}};a}_{mn}(-{\bm{\mathbf{k}}})\right)^{*}.$
(101)

##  Appendix D Computataional details

###  D.1 Details of density functional theory calculations

We used the Quantum ESPRESSO package Giannozzi  _et al._ (2017) to perform
density functional theory calculations. For the self-consistent field
calculation, we used an unshifted 16$\times$16$\times$12 $k$-point grid, a
kinetic energy cutoff of 80 Ry, fully relativistic ONCV pseudopotentials
Hamann (2013) taken from the PseudoDojo library (v0.4) van Setten  _et al._
(2018), and the Perdew-Burke-Ernzerhof (PBE) exchange-correlation functional
Perdew  _et al._ (1996). The lattice parameters and atomic coordinates were
optimized until the stresses and forces were less than $1.0\times 10^{-6}$
Ry/Bohr3 and $6.0\times 10^{-4}$ Ry/Bohr, respectively. The optimized lattice
parameters were $a$=4.44 Å and $c$=7.39 Å, consistent with those of Ref.
Brousseau-Couture  _et al._ (2020). Although the PBE functional overestimates
lattice parameters compared to the experimental values Ishizaka  _et al._
(2011), it gives band gap and Rashba splitting in reasonable agreement with
experiment Brousseau-Couture  _et al._ (2020).

We used the Wannier90 package Pizzi  _et al._ (2020) to construct the Wannier-
function-based tight-binding model. The Brillouin zone was sampled with an
8$\times$8$\times$6 grid for Wannierization. We construct 18 Wannier functions
using the atomic $p$ orbitals as initial guesses. Disentanglement was not used
since the target bands are isolated from other bands. To preserve the crystal
symmetries, we did not perform the maximal localization step. The centers of
the Wannier functions for the tight-binding model were calculated using the
translationally invariant formula: Eq. (31) of Ref. Marzari and Vanderbilt
(1997).

###  D.2 Details of photoconductivity calculations

We used a modified version of Wannier90 Pizzi  _et al._ (2020) for the
Wannier-interpolation calculation of the photoconductivity. We sampled the
Brillouin zone using a grid shifted by half the grid spacing along all axes.
Using the shifted grid speeded up the convergence. The photoconductivity was
calculated using a dense $k$-point grid of 800$\times$800$\times$800 to obtain
converged values. A fixed numerical smearing of 20 meV was applied to the
delta functions and principal values involving the energy difference between
the initial and final states. A detailed convergence study is shown in Figs.
9-11. The Fermi-Dirac occupation was calculated at zero temperature.

To avoid numerical problems related to near-degenerate states, we regularized
the denominator including intermediate states, such as $1/\omega_{np}$ in Eq.
(5), to $\real[1/(\omega_{np}+i\eta)]$ Nastos and Sipe (2006); Ibañez-Azpiroz
_et al._ (2018). The broadening parameter was set to $\eta=1~{}\mathrm{meV}$
unless noted otherwise.

We use the diagonal tight-binding approximation (TBA) within which the
position matrix elements between different Wannier functions are neglected
Ibañez-Azpiroz  _et al._ (2018). This approximation is needed because using
the existing Wannier interpolation methods, it is theoretically impossible to
calculate the spin-velocity derivative [Eq. (5)] without erroneously
truncating the summation over bands. By using the diagonal tight-binding
approximation, the sum rules are modified to include only a finite number of
bands Ibañez-Azpiroz  _et al._ (2018) and are satisfied without any truncation
error.

The accuracy of the diagonal tight-binding approximation can be tested by
calculating the charge shift current and charge and spin injection currents
without the approximation. Such a calculation is possible because the
injection current [Eq. (10)] does not involve the spin-velocity derivative,
and a truncation-error-free expression exists for the charge shift current
Ibañez-Azpiroz  _et al._ (2018). In Figs. 7 and 8, we show that the diagonal
tight-binding approximation changes the photoconductivity only slightly.

##  Appendix E Additional computational results

Figure 5:  (a-c) Spin $S^{x}$ and (d-f) spin $S^{z}$ photoconductivity of
BiTeI. The vertical dashed lines indicate the direct band gap of 0.35 eV.
Figure 6:  Spin polarization of the subgap pure spin current along the $z$
direction. The arrows indicate the in-plane spin polarization of the current.
The $S^{z}$ component is zero due to symmetry. Light frequency is
$\Omega=0.30$ eV. The scale of the arrow is smaller than in Fig. 4.  Figure 7:
Charge shift current calculated with and without the diagonal tight-binding
approximation (diag. TBA).  Figure 8:  Spin $S^{z}$ injection current
calculated with and without the diagonal tight-binding approximation.  Figure
9:  Convergence of photoconductivity with respect to the size of the $k$-point
grid. (a-c) Charge and (d-f) spin $S^{z}$ photoconductivity calculated with an
$N_{k}\times N_{k}\times N_{k}$ $k$-point grid. The plotted photoconductivity
tensor components are those plotted in Fig. 3(B-G).  Figure 10:  Convergence
of photoconductivity with respect to $\eta$. (a-c) Charge and (d-f) spin
$S^{z}$ photoconductivity calculated with $\eta=$1 meV and 0.1 meV. The
plotted photoconductivity tensor components are those plotted in Fig. 3(B-G).
Figure 11:  Convergence of charge rectification conductivity. The maximum
absolute value is calculated in the 0 eV to 2 eV window. The size of the
$k$-point grid is $N_{k}\times N_{k}\times N_{k}$. The lines are a guide to
the eye.

Figure 5 shows the spin $S^{x}$ and $S^{z}$ photoconductivity for currents
flowing along the $x$ direction.

Figure 6 shows the spin polarization of the subgap spin conductivity for
currents flowing along the $z$ direction as a function of light polarization
angles. Note that the size of the spin conductivity shown in Fig. 6 is an
order of magnitude smaller than that of Fig. 4. The spin photocurrent flowing
along the $z$ direction in BiTeI is small because the corresponding group
velocity of the carriers is small in BiTeI, a layered compound.

In Fig. 7 and Fig. 8, we show the charge shift and spin injection currents
calculated with and without the diagonal tight-binding approximation. We find
that the diagonal tight-binding approximation changes the photocurrent only
slightly.

Figures 9-11 show the convergence of the photoconductivity with respect to the
$k$-point grid size. In Fig. 9, we find that all photoconductivity elements
are already converged at a $400\times 400\times 400$ $k$-point grid. Figure 10
shows that all photoconductivity elements are converged with respect to $\eta$
at $\eta=1$ meV. In Fig. 11, we show that the charge rectification current
converges to zero in the limit of an infinitely fine $k$-point grid as
$\eta\rightarrow 0$, as expected in insulators.

## References

  * Grinberg  _et al._ (2013) Ilya Grinberg, D. Vincent West, Maria Torres, Gaoyang Gou, David M. Stein, Liyan Wu, Guannan Chen, Eric M. Gallo, Andrew R. Akbashev, Peter K. Davies, Jonathan E. Spanier, and Andrew M. Rappe, “Perovskite oxides for visible-light-absorbing ferroelectric and photovoltaic materials,” Nature 503, 509–512 (2013).
  * Bonaccorso  _et al._ (2010) Francesco Bonaccorso, Z Sun, TA Hasan, and AC Ferrari, “Graphene photonics and optoelectronics,” Nature photonics 4, 611 (2010).
  * Bhat  _et al._ (2005) R. D. R Bhat, F. Nastos, Ali Najmaie, and J. E. Sipe, “Pure Spin Current from One-Photon Absorption of Linearly Polarized Light in Noncentrosymmetric Semiconductors,” Physical Review Letters 94, 096603 (2005).
  * Ivchenko and Tarasenko (2008) E L Ivchenko and S A Tarasenko, “Pure spin photocurrents,” Semiconductor Science and Technology 23, 114007 (2008).
  * Young  _et al._ (2013) Steve M. Young, Fan Zheng, and Andrew M. Rappe, “Prediction of a Linear Spin Bulk Photovoltaic Effect in Antiferromagnets,” Physical Review Letters 110, 057201 (2013).
  * Tarasenko and Ivchenko (2005) Sergei Anatol’evich Tarasenko and Eugeniyus Levovich Ivchenko, “Pure spin photocurrents in low-dimensional structures,” Journal of Experimental and Theoretical Physics Letters 81, 231–235 (2005).
  * Zhao  _et al._ (2005) Hui Zhao, Xinyu Pan, Arthur L. Smirl, R. D. R. Bhat, Ali Najmaie, J. E. Sipe, and H. M. van Driel, “Injection of ballistic pure spin currents in semiconductors by a single-color linearly polarized beam,” Phys. Rev. B 72, 201302 (2005).
  * Ganichev  _et al._ (2006) Sergey D. Ganichev, Vasily V. Bel’kov, Sergey A. Tarasenko, Sergey N. Danilov, Stephan Giglberger, Christoph Hoffmann, Eougenious L. Ivchenko, Dieter Weiss, Werner Wegscheider, Christian Gerl, Dieter Schuh, Joachim Stahl, Jo De Boeck, Gustaaf Borghs, and Wilhelm Prettl, “Zero-bias spin separation,” Nature Physics 2, 609–613 (2006).
  * Fei  _et al._ (2020) Ruixiang Fei, Xiaobo Lu, and Li Yang, “Intrinsic Spin Photogalvanic Effect in Nonmagnetic Insulator,” arXiv:2006.10690 [cond-mat, physics:physics] (2020), arXiv:2006.10690 [cond-mat, physics:physics] .
  * Xiao  _et al._ (2021) Rui-Chun Xiao, Ding-Fu Shao, Yu-Hang Li, and Hua Jiang, “Spin photogalvanic effect in two-dimensional collinear antiferromagnets,” npj Quantum Materials 6, 1–6 (2021).
  * Ganichev  _et al._ (2001) S. D. Ganichev, E. L. Ivchenko, S. N. Danilov, J. Eroms, W. Wegscheider, D. Weiss, and W. Prettl, “Conversion of Spin into Directed Electric Current in Quantum Wells,” Physical Review Letters 86, 4358–4361 (2001).
  * Hamamoto  _et al._ (2017) Keita Hamamoto, Motohiko Ezawa, Kun Woo Kim, Takahiro Morimoto, and Naoto Nagaosa, “Nonlinear spin current generation in noncentrosymmetric spin-orbit coupled systems,” Physical Review B 95, 224430 (2017).
  * Cui  _et al._ (2007) Xiao-Dong Cui, Shun-Qing Shen, Jian Li, Yang Ji, Weikun Ge, and Fu-Chun Zhang, “Observation of electric current induced by optically injected spin current,” Applied Physics Letters 90, 242115 (2007).
  * Fei  _et al._ (2021) Ruixiang Fei, Shuaiqin Yu, Yan Lu, Linghan Zhu, and Li Yang, “Switchable Enhanced Spin Photocurrent in Rashba and Cubic Dresselhaus Ferroelectric Semiconductors,” Nano Letters 21, 2265–2271 (2021).
  * Kim  _et al._ (2017) Kun Woo Kim, Takahiro Morimoto, and Naoto Nagaosa, “Shift charge and spin photocurrents in Dirac surface states of topological insulator,” Physical Review B 95, 035134 (2017).
  * Xu  _et al._ (2021) Haowei Xu, Hua Wang, Jian Zhou, and Ju Li, “Pure spin photocurrent in non-centrosymmetric crystals: bulk spin photovoltaic effect,” Nature Communications 12, 1–9 (2021).
  * Mu  _et al._ (2021) Xingchi Mu, Yiming Pan, and Jian Zhou, “Pure bulk orbital and spin photocurrent in two-dimensional ferroelectric materials,” npj Computational Materials 7, 61 (2021).
  * Železný  _et al._ (2017) Jakub Železný, Yang Zhang, Claudia Felser, and Binghai Yan, “Spin-polarized current in noncollinear antiferromagnets,” Phys. Rev. Lett. 119, 187204 (2017).
  * Zhang  _et al._ (2018) Yang Zhang, Jakub Železný, Yan Sun, Jeroen van den Brink, and Binghai Yan, “Spin Hall effect emerging from a noncollinear magnetic lattice without spin–orbit coupling,” New Journal of Physics 20, 073028 (2018).
  * Sipe and Shkrebtii (2000) J. E. Sipe and A. I. Shkrebtii, “Second-order optical response in semiconductors,” Physical Review B 61, 5337–5352 (2000).
  * Zhang  _et al._ (2019) Yang Zhang, Tobias Holder, Hiroaki Ishizuka, Fernando de Juan, Naoto Nagaosa, Claudia Felser, and Binghai Yan, “Switchable magnetic bulk photovoltaic effect in the two-dimensional magnet CrI3,” Nature Communications 10, 3783 (2019).
  * de Juan  _et al._ (2020) F. de Juan, Y. Zhang, T. Morimoto, Y. Sun, J. E. Moore, and A. G. Grushin, “Difference frequency generation in topological semimetals,” Physical Review Research 2, 012017 (2020).
  * Holder  _et al._ (2020) Tobias Holder, Daniel Kaplan, and Binghai Yan, “Consequences of time-reversal-symmetry breaking in the light-matter interaction: Berry curvature, quantum metric, and diabatic motion,” Physical Review Research 2, 033100 (2020).
  * Ahn  _et al._ (2020) Junyeong Ahn, Guang-Yu Guo, and Naoto Nagaosa, “Low-frequency divergence and quantum geometry of the bulk photovoltaic effect in topological semimetals,” Phys. Rev. X 10, 041041 (2020).
  * Gao  _et al._ (2020) Lingyuan Gao, Zachariah Addison, E. J. Mele, and Andrew M. Rappe, “Intrinsic Fermi Surface Contribution to the Circular Photogalvanic Effect,” arXiv:2011.06542 [cond-mat] (2020), arXiv:2011.06542 [cond-mat] .
  * Watanabe and Yanase (2021) Hikaru Watanabe and Youichi Yanase, “Chiral Photocurrent in Parity-Violating Magnet and Enhanced Response in Topological Antiferromagnet,” Physical Review X 11, 011001 (2021).
  * Ventura  _et al._ (2017) G. B. Ventura, D. J. Passos, J. M. B. Lopes dos Santos, J. M. Viana Parente Lopes, and N. M. R. Peres, “Gauge covariances and nonlinear optical responses,” Physical Review B 96, 035431 (2017).
  * Parker  _et al._ (2019) Daniel E. Parker, Takahiro Morimoto, Joseph Orenstein, and Joel E. Moore, “Diagrammatic approach to nonlinear optical response with application to Weyl semimetals,” Physical Review B 99, 045121 (2019).
  * Passos  _et al._ (2018) D. J. Passos, G. B. Ventura, J. M. Viana Parente Lopes, J. M. B. Lopes dos Santos, and N. M. R. Peres, “Nonlinear optical responses of crystalline systems: Results from a velocity gauge analysis,” Physical Review B 97, 235446 (2018).
  * de Juan  _et al._ (2017) Fernando de Juan, Adolfo G. Grushin, Takahiro Morimoto, and Joel E Moore, “Quantized circular photogalvanic effect in Weyl semimetals,” Nature Communications 8, 15995 (2017).
  * Guo  _et al._ (2008) G. Y. Guo, S. Murakami, T.-W. Chen, and N. Nagaosa, “Intrinsic Spin Hall Effect in Platinum: First-Principles Calculations,” Physical Review Letters 100, 096401 (2008).
  * von Baltz and Kraut (1981) Ralph von Baltz and Wolfgang Kraut, “Theory of the bulk photovoltaic effect in pure crystals,” Phys. Rev. B 23, 5590–5596 (1981).
  * Jahn (1949) HA Jahn, “Note on the bhagavantam–suranarayana method of enumerating the physical constants of crystals,” Acta Crystallographica 2, 30–33 (1949).
  * Gallego  _et al._ (2019) S. V. Gallego, J. Etxebarria, L. Elcoro, E. S. Tasci, and J. M. Perez-Mato, “Automatic calculation of symmetry-adapted tensors in magnetic and non-magnetic materials: A new tool of the Bilbao Crystallographic Server,” Acta Cryst. A 75, 438–447 (2019).
  * Momma and Izumi (2011) Koichi Momma and Fujio Izumi, “Vesta 3 for three-dimensional visualization of crystal, volumetric and morphology data,” Journal of applied crystallography 44, 1272–1276 (2011).
  * Ishizaka  _et al._ (2011) K. Ishizaka, M. S. Bahramy, H. Murakawa, M. Sakano, T. Shimojima, T. Sonobe, K. Koizumi, S. Shin, H. Miyahara, A. Kimura, K. Miyamoto, T. Okuda, H. Namatame, M. Taniguchi, R. Arita, N. Nagaosa, K. Kobayashi, Y. Murakami, R. Kumai, Y. Kaneko, Y. Onose, and Y. Tokura, “Giant Rashba-type spin splitting in bulk BiTeI,” Nature Materials 10, 521–526 (2011).
  * Marzari and Vanderbilt (1997) Nicola Marzari and David Vanderbilt, “Maximally-localized generalized Wannier functions for composite energy bands,” Physical Review B 56, 12847–12865 (1997).
  * Saitoh  _et al._ (2006) E. Saitoh, M. Ueda, H. Miyajima, and G. Tatara, “Conversion of spin current into charge current at room temperature: Inverse spin-hall effect,” Applied Physics Letters 88, 182509 (2006), https://doi.org/10.1063/1.2199473 .
  * Tan  _et al._ (2016) Liang Z Tan, Fan Zheng, Steve M Young, Fenggong Wang, Shi Liu, and Andrew M Rappe, “Shift current bulk photovoltaic effect in polar materials—hybrid and oxide perovskites and beyond,” npj Computational Materials 2, 16026 (2016).
  * Cook  _et al._ (2017) Ashley M. Cook, Benjamin M. Fregoso, Fernando de Juan, Sinisa Coh, and Joel E. Moore, “Design principles for shift current photovoltaics,” Nature Communications 8, 14176 (2017).
  * Rangel  _et al._ (2017) Tonatiuh Rangel, Benjamin M. Fregoso, Bernardo S. Mendoza, Takahiro Morimoto, Joel E. Moore, and Jeffrey B. Neaton, “Large Bulk Photovoltaic Effect and Spontaneous Polarization of Single-Layer Monochalcogenides,” Physical Review Letters 119, 067402 (2017).
  * Ivchenko  _et al._ (1990) EL Ivchenko, Yu B Lyanda-Geller, and GE Pikus, “Current of thermalized spin-oriented photocarriers,” Sov. Phys.-JETP 71, 550–557 (1990).
  * Ganichev  _et al._ (2002) SD Ganichev, EL Ivchenko, VV Bel’Kov, SA Tarasenko, M Sollinger, Dieter Weiss, Werner Wegscheider, and Wilhelm Prettl, “Spin-galvanic effect,” Nature 417, 153–156 (2002).
  * Budkin and Tarasenko (2020) G V Budkin and S A Tarasenko, “Thermal generation of shift electric current,” New Journal of Physics 22, 013005 (2020).
  * Kaplan  _et al._ (2020) Daniel Kaplan, Tobias Holder, and Binghai Yan, “Nonvanishing Subgap Photocurrent as a Probe of Lifetime Effects,” Physical Review Letters 125, 227401 (2020).
  * Lihm and Park (2021) Jae-Mo Lihm and Cheol-Hwan Park, “Wannier function perturbation theory: Localized representation and interpolation of wave function perturbation,” Phys. Rev. X 11, 041053 (2021).
  * Shi  _et al._ (2006) Junren Shi, Ping Zhang, Di Xiao, and Qian Niu, “Proper definition of spin current in spin-orbit coupled systems,” Phys. Rev. Lett. 96, 076604 (2006).
  * Sun  _et al._ (2008) Qing-feng Sun, X. C. Xie, and Jian Wang, “Persistent spin current in nanodevices and definition of the spin current,” Phys. Rev. B 77, 035327 (2008).
  * Qiao  _et al._ (2018) Junfeng Qiao, Jiaqi Zhou, Zhe Yuan, and Weisheng Zhao, “Calculation of intrinsic spin Hall conductivity by Wannier interpolation,” Physical Review B 98, 214402 (2018).
  * Wang  _et al._ (2020) Haihang Wang, Priya Gopal, Silvia Picozzi, Stefano Curtarolo, Marco Buongiorno Nardelli, and Jagoda Sławińska, “Spin Hall effect in prototype Rashba ferroelectrics GeTe and SnTe,” npj Computational Materials 6, 7 (2020).
  * Rashba (2003) Emmanuel I. Rashba, “Spin currents in thermodynamic equilibrium: The challenge of discerning transport currents,” Phys. Rev. B 68, 241315 (2003).
  * Blount (1962) E.I. Blount, “Formalisms of band theory,” (Academic Press, 1962) pp. 305 – 373.
  * Giannozzi  _et al._ (2017) P Giannozzi, O Andreussi, T Brumme, O Bunau, M Buongiorno Nardelli, M Calandra, R Car, C Cavazzoni, D Ceresoli, M Cococcioni, N Colonna, I Carnimeo, A Dal Corso, S de Gironcoli, P Delugas, R A DiStasio, A Ferretti, A Floris, G Fratesi, G Fugallo, R Gebauer, U Gerstmann, F Giustino, T Gorni, J Jia, M Kawamura, H-Y Ko, A Kokalj, E Küçükbenli, M Lazzeri, M Marsili, N Marzari, F Mauri, N L Nguyen, H-V Nguyen, A Otero-de-la-Roza, L Paulatto, S Poncé, D Rocca, R Sabatini, B Santra, M Schlipf, A P Seitsonen, A Smogunov, I Timrov, T Thonhauser, P Umari, N Vast, X Wu, and S Baroni, “Advanced capabilities for materials modelling with Quantum ESPRESSO,” Journal of Physics: Condensed Matter 29, 465901 (2017).
  * Hamann (2013) D. R. Hamann, “Optimized norm-conserving Vanderbilt pseudopotentials,” Physical Review B 88, 085117 (2013).
  * van Setten  _et al._ (2018) M.J. van Setten, M. Giantomassi, E. Bousquet, M.J. Verstraete, D.R. Hamann, X. Gonze, and G.-M. Rignanese, “The PseudoDojo: Training and grading a 85 element optimized norm-conserving pseudopotential table,” Computer Physics Communications 226, 39–54 (2018).
  * Perdew  _et al._ (1996) John P. Perdew, Kieron Burke, and Matthias Ernzerhof, “Generalized gradient approximation made simple,” Phys. Rev. Lett. 77, 3865–3868 (1996).
  * Brousseau-Couture  _et al._ (2020) Veronique Brousseau-Couture, Gabriel Antonius, and Michel Cote, “Temperature dependence of the topological phase transition of BiTeI from first principles,” Physical Review Research 2, 023185 (2020).
  * Pizzi  _et al._ (2020) Giovanni Pizzi, Valerio Vitale, Ryotaro Arita, Stefan Blügel, Frank Freimuth, Guillaume Géranton, Marco Gibertini, Dominik Gresch, Charles Johnson, Takashi Koretsune, Julen Ibañez-Azpiroz, Hyungjun Lee, Jae-Mo Lihm, Daniel Marchand, Antimo Marrazzo, Yuriy Mokrousov, Jamal I Mustafa, Yoshiro Nohara, Yusuke Nomura, Lorenzo Paulatto, Samuel Poncé, Thomas Ponweiser, Junfeng Qiao, Florian Thöle, Stepan S Tsirkin, Małgorzata Wierzbowska, Nicola Marzari, David Vanderbilt, Ivo Souza, Arash A Mostofi, and Jonathan R Yates, “Wannier90 as a community code: New features and applications,” Journal of Physics: Condensed Matter 32, 165902 (2020).
  * Nastos and Sipe (2006) F. Nastos and J. E. Sipe, “Optical rectification and shift currents in gaas and gap response: Below and above the band gap,” Phys. Rev. B 74, 035201 (2006).
  * Ibañez-Azpiroz  _et al._ (2018) Julen Ibañez-Azpiroz, Stepan S. Tsirkin, and Ivo Souza, “ _Ab Initio_ calculation of the shift photocurrent by Wannier interpolation,” Physical Review B 97, 245143 (2018).