# Evolution of $\Lambda$ polarization in the hadronic phase of heavy-ion
collisions

Yifeng Sun  Department of Physics and Astronomy, University of Catania, Via S.
Sofia 64, 1-95125 Catania, Italy  Laboratori Nazionali del Sud, INFN-LNS, Via
S. Sofia 62, I-95123 Catania, Italy  Zhen Zhang<EMAIL_ADDRESS>Sino-French Institute of Nuclear Engineering and Technology, Sun Yat-sen
University, Zhuhai 519082, China  Che Ming Ko  Cyclotron Institute and
Department of Physics and Astronomy, Texas A$\&$M University, College Station,
Texas 77843-3366, USA  Wenbin Zhao  Department of Physics and Astronomy, Wayne
State University, Detroit, Michigan 48201, USA

###### Abstract

Using the AMPT + MUSIC+UrQMD hybrid model, we study the global and local spin
polarizations of $\Lambda$ hyperons as functions of the freeze-out temperature
of the spin degree of freedom in the hadronic phase of Au+Au collisions at
$\sqrt{s_{NN}}=19.6$ GeV. Including contributions from both the thermal
vorticity and thermal shear of the hadronic matter, we find that, with the
spin freeze-out temperature dropping from the hadronization temperature of 160
MeV to 110 MeV at the kinetic freeze-out, both the global and local spin
polarizations of $\Lambda$ hyperons due to the thermal vorticity decrease by a
factor of two, while those due to the thermal shear decrease quickly and
become negligibly small at 140 MeV. Our results suggest the importance of
understanding the dynamical evolution of the spin degree of freedom in the
hadronic stage in relativistic heavy-ion collisions.

##  I Introduction

In non-central heavy ion collisions at relativistic energies, a large amount
of orbital angular momentum from the two colliding nuclei are transferred to
the produced quark-gluon plasma (QGP), creating thus the most vortical fluid
of the order of $10^{21}-10^{22}$ s-1 in known physical systems Adamczyk  _et
al._ (2017); Deng and Huang (2016); Jiang  _et al._ (2016); Csernai  _et al._
(2013, 2014); Baznat  _et al._ (2013); Teryaev and Usubov (2015). Due to their
spin-orbit interactions, quarks in the QGP and $\Lambda$ hyperons formed after
the hadronization become polarized along the direction of the total orbital
angular momentum Liang and Wang (2006); Gao  _et al._ (2008), as observed in
experiments by the STAR Collaboration Adamczyk  _et al._ (2017); Adam  _et
al._ (2021). Because of the non-uniformity of the vorticity field, local
structures in the $\Lambda$ spin polarization have also been found both in
theoretical studies Pang  _et al._ (2016); Becattini and Karpenko (2018); Xia
_et al._ (2018) and in experimental measurements Adam  _et al._ (2019).
Although various theoretical models, such as those based on the hydrodynamic
approach Becattini  _et al._ (2013); Karpenko and Becattini (2017); Xie  _et
al._ (2016, 2017); Becattini  _et al._ (2015); Wu  _et al._ (2019); Fu  _et
al._ (2021a), the transport approach Li  _et al._ (2017); Kolomeitsev  _et
al._ (2018); Xia  _et al._ (2018); Wei  _et al._ (2019); Deng  _et al._ (2021)
and the non-equilibrium chiral kinetic approach Sun and Ko (2017, 2019); Liu
_et al._ (2020), have successfully described the measured $\Lambda$ global
spin polarization, most of them have failed to explain the measured $\Lambda$
local spin polarizations. A plausible explanation of the latter has been
provided in the chiral kinetic approach through the induced quadrupole axial
charge distribution in the transverse plane of a heavy ion collision Sun and
Ko (2019); Liu  _et al._ (2020). Also, it has recently been pointed out that
adding the thermal shear contribution to that from the thermal vorticity in
the fluid dynamic approach can potentially describe the measured azimuthal
angle dependence of $\Lambda$ local spin polarizations Liu and Yin (2021);
Becattini  _et al._ (2021a); Yi  _et al._ (2021); Liu and Huang (2021).
Indeed, the correct azimuthal angle dependence is obtained in the “strange
memory” scenario of Ref. Fu  _et al._ (2021b), which assumes that the
$\Lambda$ spin polarization is identical to its strange quark spin
polarization like in the quark coalescence model for $\Lambda$ production Sun
and Ko (2017), after including contributions from both thermal vorticity and
thermal shear. A similar result can also be achieved in the “isothermal local
equilibrium” scenario of Ref. Becattini  _et al._ (2021b), in which the
hadronization is assumed to take place at a constant temperature to eliminate
the contribution from its space-time gradients in the calculation of the
thermal vorticity and shear Becattini  _et al._ (2021b).

In almost all these theoretical studies, the $\Lambda$ spin polarization in a
heavy ion collision is calculated at the end of the partonic phase and
compared to experimental measurements. Since the strength of the vorticity
field decreases as the hadronic matter expands Jiang  _et al._ (2016), it is
of interest to study how the $\Lambda$ polarization changes during the
hadronic evolution of heavy ion collisions. In this work, we study the time
evolution of the global and local $\Lambda$ spin polarizations during the
expansion of the hadronic matter produced in non-central Au+Au collisions at
$\sqrt{s_{NN}}=19.6$ GeV by using the MUSIC+UrQMD hybrid model. Because of the
lack of a dynamical description of the spin degree of freedom in the hadronic
transport model, we study the dependence of $\Lambda$ spin polarizations on
the temperature by assuming that they are in thermal equilibrium with the
vorticity field in the expanding hadronic matter.

The paper is organized as follows: In Sec. II, we discuss the relations
between the spin polarization of a spin-1/2 fermion and the thermal vorticity
and shear in a hadronic matter. We then present the numerical results from our
study on the global and local spin polarization of $\Lambda$ hyperons at
different freeze-out temperatures for the spin degree of freedom. Finally, a
brief conclusion and discussion is given in Sec. IV.

##  II Spin polarization of a fermion in thermal equilibrium

The spin polarization of a fermion in a thermal medium of temperature $T$
depends on the four-temperature vector $\beta=u/T$, where $u$ is the four flow
velocity of the local medium. To the leading order in the gradient of $\beta$,
the spin polarization vector $S$ of a fermion of four momentum $p$ at four
space-time coordinate $x$ can be written as Liu and Yin (2021); Becattini  _et
al._ (2021a); Yi  _et al._ (2021); Liu and Huang (2021)

$\displaystyle S^{\mu}(x,p)=$
$\displaystyle-\frac{1}{8m}(1-n_{F})\epsilon^{\mu\nu\rho\sigma}p_{v}\varpi_{\rho\sigma}(x)$
(1)
$\displaystyle-\frac{1}{4m}(1-n_{F})\epsilon^{\mu\nu\rho\sigma}p_{v}\frac{n_{\rho}p^{\lambda}\xi_{\lambda\sigma}(x)}{n\cdot
p}.$

In the above,

$\displaystyle\varpi_{\rho\sigma}=\frac{1}{2}(\partial_{\sigma}\beta_{\rho}-\partial_{\rho}\beta_{\sigma}),\quad\xi_{\rho\sigma}=\frac{1}{2}(\partial_{\sigma}\beta_{\rho}+\partial_{\rho}\beta_{\sigma})$
(2)

are the thermal vorticity and thermal shear of the medium, respectively,
$n_{F}$ is the Fermi-Dirac distribution function, and $n$ is a unit four
vector that specifies the frame of reference. Although different forms of
$S^{\mu}$ are used in Refs. Liu and Yin (2021); Becattini  _et al._ (2021a);
Yi  _et al._ (2021); Liu and Huang (2021), they all become the same if $n$ is
taken to be the flow field $u$ in the fluid, which we adopt in this work. Our
calculation using Eq. (2) is similar to that in Ref. Fu  _et al._ (2021b) but
different from that based on the assumption of isothermal local equilibrium
adopted in Ref. Becattini  _et al._ (2021b), where both the spatial and
temporal gradients of the temperature are neglected, i.e., only contributions
from the kinetic vorticity and shear are included.

By decomposing the two terms in Eq.(1) into the following two components,

$\displaystyle\bm{\varpi}_{T}=\frac{1}{2}\left[\bm{\nabla}\left(\frac{u^{0}}{T}\right)+\partial_{t}\left(\frac{\bm{u}}{T}\right)\right]-u^{0}\bm{f}+\bm{u}f^{0},$
$\displaystyle\bm{\varpi}_{S}=\frac{1}{2}\bm{\nabla}\times\left(\frac{\bm{u}}{T}\right)+\bm{u}\times\bm{f},$
(3)

with

$\displaystyle u^{0}=\gamma,\quad\bm{u}=\gamma\bm{v},$ $\displaystyle
f^{0}=\frac{p^{\lambda}\xi_{0\lambda}}{u\cdot
p},\quad\bm{f}=-(f_{1},f_{2},f_{3}),\quad
f_{i}=\frac{p^{\lambda}\xi_{i\lambda}}{u\cdot p},$ (4)

Eq. (1) can be rewritten as

$\displaystyle S^{0}(x,p)=\frac{1}{4m}\bm{p}\cdot\bm{\varpi}_{S},$
$\displaystyle\bm{S}(x,p)=\frac{1}{4m}(E_{p}\bm{\varpi}_{S}+\bm{p}\times\bm{\varpi}_{T}),$
(5)

where $E_{p}=\sqrt{m^{2}+{\bf p}^{2}}$ and ${\bf p}$ are, respectively, the
time and space components of the four momentum $p$.

Boosting the spin vector $S$ in the fluid frame to the particle’s rest frame
and using the fact that the spin of a fermion equal to $1/2$, one finds the
spin polarization ${\bf P}$ measured in experiments to be

$\displaystyle\bm{P}=2\bm{S}^{*}=2\left[\bm{S}-\frac{\bm{p}\times\bm{S}}{E_{p}(E_{p}+m)}\bm{p}\right].$
(6)

##  III Numerical results

To study the evolution of $\Lambda$ spin polarization in the hadronic phase of
a heavy ion collision, we first use the MUSIC hydrodynamic model Schenke  _et
al._ (2011) with the initial conditions taken from A Multi-Phase Transport
(AMPT) model Lin  _et al._ (2005) to simulate the QGP phase of the collision.
After converting the fluid elements on the freeze-out hypersurface into
hadrons by the Cooper-Frye formula, we adopt the UrQMD model Bass  _et al._
(1998); Bleicher  _et al._ (1999) to simulate the evolution of these hadrons
due to their scatterings. This hybrid approach has successfully described
various soft hadronic observables, such as the charged particle yields as well
as their transverse momentum spectra and flow anisotropies in relativistic
heavy ion collisions at both RHIC and LHC energies Shen and Alzhrani (2020);
Zhao  _et al._ (2020); Schenke  _et al._ (2020).

To determine the local temperature and flow field as well as their gradients
in time and space in the hadronic matter, we use the coarse-grained method
with the time step $\Delta t=0.5$ fm$/c$, and the cell size $\Delta x=\Delta
y=0.5$ fm and $\Delta\eta=0.2$, where $\eta=\frac{1}{2}\ln\frac{t+z}{t-z}$ is
the space-time rapidity. In this method, the temperature $T$ of a local cell
is calculated from the energy density in the cell at its rest frame by using
the equation of state from the lattice QCD calculations Borsanyi  _et al._
(2014). For the flow field $\bm{v}$ in each cell, it is calculated from the
average velocity of the $N$ hadrons in the cell according to

$\displaystyle\bm{v}(t,x,y,z)=\frac{\sum_{i}\frac{\bm{p_{i}}}{E_{i}}}{N},$ (7)

where ${\bf p}_{i}$ and $E_{i}$ are the momentum and energy of the $i$-th
hadron in the cell.

Figure 1: (Color online) Spin freeze-out temperature dependence of $\Lambda$
hyperon spin polarization along the total angular momentum direction due to
the thermal vorticity and the thermal shear in an expanding hadronic matter.

To illustrate how $\Lambda$ spin polarizations are affected by the temperature
at which the spin degrees of freedom freeze out in relativistic heavy ion
collisions, we consider in the present study Au+Au collisions at
$\sqrt{s_{NN}}=19.6$ GeV and 20-60% centrality. As in the experimental
measurements at RHIC by the STAR Collaboration Adam  _et al._ (2019), we only
include $\Lambda$ hyperons of $|\eta|<1$ and $p_{T}>$0.15 GeV$/c$ in our
analysis. In Fig. 1, we first show by the solid red line the global spin
polarization of $\Lambda$ ($\overline{\Lambda}$) hyperons calculated from the
thermal vorticity along the direction of the total angular momentum of the
hadronic matter as a function of the local temperature in the hadronic matter.
It is seen that if the spin degree of freedom freezes out at $T_{s}=160$ MeV,
which is usually assumed in the literature Li  _et al._ (2017); Sun and Ko
(2017); Fu  _et al._ (2021a) as mentioned in the introduction, the $\Lambda$
global polarization is 5.3$\times 10^{-3}$ and is close to the value in Ref.
Fu  _et al._ (2021a) based on the same AMPT+MUSIC hybrid approach without the
UrQMD afterburner. With decreasing temperature as the hadronic matter expands,
the $\Lambda$ global spin polarization decreases to 2.5$\times 10^{-3}$ at the
kinetic freeze-out temperature of $T_{s}=110$ MeV, when hadron momentum
spectra and anisotropies stop changing. The above result thus shows that the
$\Lambda$ global spin polarization has a strong temperature dependence in the
expanding hadronic matter if it continues to be in thermal equilibrium .

Shown by the blue dashed line in Fig. 1 is the $\Lambda$ global spin
polarization generated by the thermal shear in the hadronic matter, i.e., the
second term in the RHS of Eq. (1), as a function of local temperature. It is
seen that at $T_{s}=160$ MeV, the $\Lambda$ global polarization due to the
thermal shear is about 2.9$\times 10^{-3}$, which is comparable to that due to
the thermal vorticity. We note that the total spin polarization of $\Lambda$
hyperons due to the thermal shear should have been zero without the cut of
$p_{T}\geq 0.15$ GeV/$c$ and $|\eta|<$1 in the $\Lambda$ momentum. With
decreasing spin freeze-out temperature, the $\Lambda$ global spin polarization
due to the thermal shear decreases faster compared to that due to the thermal
vorticity, and it becomes negligible at $T_{s}=140$ MeV.

Figure 2: (Color online) The azimuthal angle dependence of the local spin
polarizations $P_{x},-P_{y}$ and $P_{z}$ of $\Lambda$ hyperons generated by
thermal vorticity and thermal shear at the spin freeze-out temperatures
$T_{s}=$160, 140, and 120 MeV.

As to the local spin polarization of $\Lambda$ hyperons as a function of its
azimuthal angle in the transverse plane of a heavy ion collision, it has
recently been extensively studied Becattini and Karpenko (2018); Sun and Ko
(2019); Liu  _et al._ (2020); Adam  _et al._ (2019). It is found in Refs. Fu
_et al._ (2021b); Becattini  _et al._ (2021b) that both thermal vorticity and
thermal shear are important in determining the local spin polarization. In the
left panels (a), (b) and (c) of Fig. 2, we show, respectively, the spin
polarization $P_{x},-P_{y}$ and $P_{z}$ of $\Lambda$ hyperons generated by the
thermal vorticity (red dashed lines) and the thermal shear (blue dashed-dotted
lines) as functions of the azimuthal angle $\phi_{p}$ of the $\Lambda$
transverse momentum if the spin degree of freedom freezes out at $T_{s}=160$
MeV. It is seen that $-P_{y}$ (left-middle panel) and $P_{z}$ (left-lower
panel) have the form of $-\cos(2\phi_{p})$ and $-\sin(2\phi_{p})$,
respectively, due to the thermal vorticity and $\cos(2\phi_{p})$ and
$\sin(2\phi_{p})$, respectively, due to the thermal shear, which are similar
to the findings in Refs. Fu  _et al._ (2021b); Becattini  _et al._ (2021b); Yi
_et al._ (2021). However, such an oscillatory azimuthal angle dependence
essentially disappears in the total $\Lambda$ spin polarization after adding
the two contributions, which disagrees with the experimental measurements Adam
_et al._ (2019), unless one adopts either the “strange memory” scenario as in
Ref. Fu  _et al._ (2021b) or the “isothermal local equilibrium” scenario as in
Ref. Becattini  _et al._ (2021b). In the present work, we focus on the effect
of hadronic evolution on the $\Lambda$ spin polarization and postpone the
study of the above two scenarios for future study.

In the left-upper panel (a) of Fig. 2, we also show for the first time the
azimuthal angle dependence of $P_{x}$ generated by the thermal vorticity and
thermal shear. It is seen that $P_{x}$ has the form of $-\sin(2\phi_{p})$ due
to the thermal vorticity and the form of $\sin(2\phi_{p})$ due to the thermal
shear. The total $P_{x}$ after adding the two oscillatory contributions is,
however, negligible.

Figure 3: (Color online) Azimuthal angle dependence of $\Lambda$ spin
polarization generated by components of the thermal shear at the spin freeze-
out temperature $T_{s}=160$ MeV.

In Fig. 3, we further show the $\Lambda$ spin polarization generated by
individual components of the thermal shear. We find that the main contribution
to $-P_{y}$ is from the $\xi_{13}$ component and that to $P_{z}$ is from the
$\xi_{11}$ and $\xi_{22}$ components. According to Eqs.(3), (5) and (6),
$\xi_{13}$ contributes to $-P_{y}$ through a term of the form
$p_{x}^{2}\xi_{13}\propto(1+\cos(2\phi_{p}))\xi_{13}$. For $P_{z}$, the
contribution from $\xi_{11}$ and $\xi_{22}$ is through a term of the form
$p_{x}p_{y}(\xi_{11}-\xi_{22})\propto\sin(2\phi_{p})(\xi_{11}-\xi_{22})$,
where $(\xi_{11}-\xi_{22})$ is positive because of the positive elliptic flow
of the hadronic matter in non-central Au+Au collisions at $\sqrt{s_{NN}}=19.6$
GeV.

We have also studied the azimuthal angle dependence of the spin polarization
of $\Lambda$ hyperons generated by the thermal vorticity and thermal shear at
a lower temperature of $T_{s}=140$ MeV during the later stage of the hadronic
evolution by assuming that the spin degree of freedom remains in thermal
equilibrium. As shown by the red dashed lines in the middle panels (d), (e)
and (f) of Fig. 2, $P_{x},-P_{y}$ and $P_{z}$ generated by the thermal
vorticity at this temperature have the same azimuthal angle dependence as
those at $T_{s}=160$ MeV. Although $P_{x}$ and $-P_{y}$ become smaller at the
lower temperature, $P_{z}$ does not change much with temperature. For the
contribution from the thermal shear at $T_{s}=140$ MeV, $P_{x}$ and $-P_{y}$
become negligible, while $P_{z}$ decreases by a factor of 3 as shown by the
blue dashed-dotted lines. By adding contributions from both thermal vorticity
and thermal shear, $P_{x}$ and $P_{z}$ have the form of $-\sin(2\phi_{p})$ and
$P_{y}$ has the form of $-\cos(2\phi_{p})$. These results suggest that the
$\Lambda$ local spin polarization depends strongly on its evolution during the
hadronic phase of heavy ion collisions.

Finally, we consider an even lower temperature of $T_{s}=120$ MeV,
corresponding approximately to the kinetic freeze-out of the expanding
hadronic matter, and the results are shown in the right panels (g), (h) and
(k) of Fig. 2. It is seen from the red dashed lines that the $\Lambda$ local
spin polarizations $P_{x},-P_{y}$ and $P_{z}$ generated by the thermal
vorticity become even smaller, although they have the same $\phi_{p}$
dependence as that at higher temperatures. However, for the $\Lambda$ spin
polarization generated by the thermal shear, its values along all $x$, $y$ or
$z$ directions become negligibly small. Thus, the total $P_{x}$, $-P_{y}$ and
$P_{z}$ have the same $\phi_{p}$ dependence as that at $T_{s}=140$ MeV, albeit
somewhat smaller.

##  IV Conclusion and Discussion

Using the AMPT model initial conditions to the MUSIC hydrodynamic model for
the QGP phase, which is followed by the UrQMD model for the hadronic phase in
relativistic heavy ion collisions, we have studied in this paper the global
and local spin polarizations of $\Lambda$ hyperons as functions of the freeze-
out temperature of the spin degree of freedom in the hadronic phase. We have
found that both the $\Lambda$ global and local spin polarizations due to the
thermal vorticity decrease by a factor of 2 if the spin freeze-out temperature
drops from 160 to 140 MeV, while those due to the thermal shear already
becomes negligibly small at temperature equal to 140 MeV. This result suggests
the importance of including in theoretical studies the dynamical evolution of
the spin degrees of freedom and its freeze-out in the expanding hadronic
matter. Such a study requires a hadronic transport model that takes into
account explicitly the spin degrees of freedom of hadrons, such as in Refs. Xu
and Li (2013); Xia  _et al._ (2014, 2016); Weickgenannt  _et al._ (2021);
Sheng  _et al._ (2021) on the nucleon spin transport in heavy ion reactions at
lower energies, which we plan to study in the future.

## ACKNOWLEDGEMENTS

We thank Huichao Song for helpful discussions and suggestions. This work was
supported by INFN-SIM national project and linea di intervento 2 for HQCDyn at
DFA-Unict (Y.S.), the National Natural Science Foundation of China under Grant
No. 11905302 (Z.Z.), the U.S. Department of Energy under Contract No. DE-
SC0015266 and the Welch Foundation under Grant No. A-1358 (C.M.K.), and the
U.S. National Science Foundation (NSF) under grant numbers ACI-2004571 (W.Z.).

## References

  * Adamczyk  _et al._ (2017) L. Adamczyk _et al._ (STAR), Nature 548, 62 (2017).
  * Deng and Huang (2016) W.-T. Deng and X.-G. Huang, Phys. Rev. C 93, 064907 (2016).
  * Jiang  _et al._ (2016) Y. Jiang, Z.-W. Lin, and J. Liao, Phys. Rev. C 94, 044910 (2016).
  * Csernai  _et al._ (2013) L. P. Csernai, V. K. Magas, and D. J. Wang, Phys. Rev. C 87, 034906 (2013).
  * Csernai  _et al._ (2014) L. P. Csernai, D. J. Wang, M. Bleicher, and H. Stöcker, Phys. Rev. C 90, 021904 (2014).
  * Baznat  _et al._ (2013) M. Baznat, K. Gudima, A. Sorin, and O. Teryaev, Phys. Rev. C 88, 061901 (2013).
  * Teryaev and Usubov (2015) O. Teryaev and R. Usubov, Phys. Rev. C 92, 014906 (2015).
  * Liang and Wang (2006) Z.-T. Liang and X.-N. Wang, Phys. Rev. Lett. 96, 039901 (2006).
  * Gao  _et al._ (2008) J.-H. Gao, S.-W. Chen, W.-T. Deng, Z.-T. Liang, Q. Wang, and X.-N. Wang, Phys. Rev. C 77, 044902 (2008).
  * Adam  _et al._ (2021) J. Adam _et al._ (STAR Collaboration), Phys. Rev. Lett. 126, 162301 (2021).
  * Pang  _et al._ (2016) L.-G. Pang, H. Petersen, Q. Wang, and X.-N. Wang, Phys. Rev. Lett. 117, 192301 (2016).
  * Becattini and Karpenko (2018) F. Becattini and I. Karpenko, Phys. Rev. Lett. 120, 012302 (2018).
  * Xia  _et al._ (2018) X.-L. Xia, H. Li, Z. Tang, and Q. Wang, Phys. Rev. C 98, 024905 (2018).
  * Adam  _et al._ (2019) J. Adam _et al._ (STAR Collaboration), Phys. Rev. Lett. 123, 132301 (2019).
  * Becattini  _et al._ (2013) F. Becattini, L. P. Csernai, and D. J. Wang, Phys. Rev. C 88, 034905 (2013).
  * Karpenko and Becattini (2017) I. Karpenko and F. Becattini, Eur. Phys. J. C 77, 213 (2017).
  * Xie  _et al._ (2016) Y. L. Xie, M. Bleicher, H. Stöcker, D. J. Wang, and L. P. Csernai, Phys. Rev. C 94, 054907 (2016).
  * Xie  _et al._ (2017) Y. Xie, D. Wang, and L. P. Csernai, Phys. Rev. C 95, 031901 (2017).
  * Becattini  _et al._ (2015) F. Becattini, G. Inghirami, V. Rolando, A. Beraudo, L. Del Zanna, A. De Pace, M. Nardi, G. Pagliara, and V. Chandra, Eur. Phys. J. C 75, 406 (2015).
  * Wu  _et al._ (2019) H.-Z. Wu, L.-G. Pang, X.-G. Huang, and Q. Wang, Phys. Rev. Research 1, 033058 (2019).
  * Fu  _et al._ (2021a) B. Fu, K. Xu, X.-G. Huang, and H. Song, Phys. Rev. C 103, 024903 (2021a).
  * Li  _et al._ (2017) H. Li, L.-G. Pang, Q. Wang, and X.-L. Xia, Phys. Rev. C 96, 054908 (2017).
  * Kolomeitsev  _et al._ (2018) E. E. Kolomeitsev, V. D. Toneev, and V. Voronyuk, Phys. Rev. C 97, 064902 (2018).
  * Wei  _et al._ (2019) D.-X. Wei, W.-T. Deng, and X.-G. Huang, Phys. Rev. C 99, 014905 (2019).
  * Deng  _et al._ (2021) X.-G. Deng, X.-G. Huang, and Y.-G. Ma, (2021), arXiv:2109.09956 [nucl-th] .
  * Sun and Ko (2017) Y. Sun and C. M. Ko, Phys. Rev. C 96, 024906 (2017).
  * Sun and Ko (2019) Y. Sun and C. M. Ko, Phys. Rev. C 99, 011903 (2019).
  * Liu  _et al._ (2020) S. Y. F. Liu, Y. Sun, and C. M. Ko, Phys. Rev. Lett. 125, 062301 (2020).
  * Liu and Yin (2021) S. Y. F. Liu and Y. Yin, JHEP 07, 188 (2021).
  * Becattini  _et al._ (2021a) F. Becattini, M. Buzzegoli, and A. Palermo, Phys. Lett. B 820, 136519 (2021a), arXiv:2103.10917 [nucl-th] .
  * Yi  _et al._ (2021) C. Yi, S. Pu, and D.-L. Yang, (2021), arXiv:2106.00238 [hep-ph] .
  * Liu and Huang (2021) Y.-C. Liu and X.-G. Huang, (2021), arXiv:2109.15301 [nucl-th] .
  * Fu  _et al._ (2021b) B. Fu, S. Y. F. Liu, L. Pang, H. Song, and Y. Yin, Phys. Rev. Lett. 127, 142301 (2021b).
  * Becattini  _et al._ (2021b) F. Becattini, M. Buzzegoli, A. Palermo, G. Inghirami, and I. Karpenko, (2021b), arXiv:2103.14621 [nucl-th] .
  * Schenke  _et al._ (2011) B. Schenke, S. Jeon, and C. Gale, Phys. Rev. Lett. 106, 042301 (2011).
  * Lin  _et al._ (2005) Z.-W. Lin, C. M. Ko, B.-A. Li, B. Zhang, and S. Pal, Phys. Rev. C 72, 064901 (2005).
  * Bass  _et al._ (1998) S. A. Bass _et al._ , Prog. Part. Nucl. Phys. 41, 255 (1998).
  * Bleicher  _et al._ (1999) M. Bleicher _et al._ , J. Phys. G 25, 1859 (1999).
  * Shen and Alzhrani (2020) C. Shen and S. Alzhrani, Phys. Rev. C 102, 014909 (2020), arXiv:2003.05852 [nucl-th] .
  * Zhao  _et al._ (2020) W. Zhao, C. Shen, C. M. Ko, Q. Liu, and H. Song, Phys. Rev. C 102, 044912 (2020), arXiv:2009.06959 [nucl-th] .
  * Schenke  _et al._ (2020) B. Schenke, C. Shen, and P. Tribedy, Phys. Rev. C 102, 044905 (2020), arXiv:2005.14682 [nucl-th] .
  * Borsanyi  _et al._ (2014) S. Borsanyi, Z. Fodor, C. Hoelbling, S. D. Katz, S. Krieg, and K. K. Szabo, Phys. Lett. B 730, 99 (2014).
  * Xu and Li (2013) J. Xu and B.-A. Li, Phys. Lett. B 724, 346 (2013).
  * Xia  _et al._ (2014) Y. Xia, J. Xu, B.-A. Li, and W.-Q. Shen, Phys. Rev. C 89, 064606 (2014).
  * Xia  _et al._ (2016) Y. Xia, J. Xu, B.-A. Li, and W.-Q. Shen, Nucl. Phys. A 955, 41 (2016).
  * Weickgenannt  _et al._ (2021) N. Weickgenannt, E. Speranza, X.-l. Sheng, Q. Wang, and D. H. Rischke, Phys. Rev. D 104, 016022 (2021).
  * Sheng  _et al._ (2021) X.-L. Sheng, N. Weickgenannt, E. Speranza, D. H. Rischke, and Q. Wang, Phys. Rev. D 104, 016029 (2021).


Supersymmetry anomalies and the Wess-Zumino Model in a supergravity background

Giorgos Eleftheriou and Peter West

Department of Mathematics

King’s College, London WC2R 2LS, UK

Abstract

We briefly recall the procedure for computing the Ward Identities in the
presence of a regulator which violates the symmetry being considered. We
compute the first non-trivial correction to the supersymmetry Ward identity of
the Wess-Zumino model in the presence of background supergravity using
dimensional regularisation. We find that the result can be removed using a
finite local counter term and so there is no supersymmetry anomaly.

Supersymmetry was discovered in 1972 and, after a slow start, it has become a
dominate theme in fundamental physics. From time to time it has been suggested
that supersymmetry is violated by quantum corrections, that is, there is a
supersymmetry anomaly. One of the problems with calculations in supersymmetric
theories is that there is no regulator which both preserves supersymmetry and
gauge symmetry and this includes the most popular regulator, namely
dimensional regularisation [1,2]. Dimensional regulation preserves gauge
symmetry and as a result it played a key role in establishing the consistency
of the standard model, however, it does not preserve supersymmetry. In the
past the claims that supersymmetry has an anomaly have been resolved by
further explicit calculation of the offending Feynman diagrams which were
causing the problem, see for example [3]. However the claims have also been
disproved by establishing that the supersymmetry Ward identities really do
hold if one took account of the fact that the regulator breaks supersymmetry
[4,5]. Indeed in this way it was shown that supersymmetry was not broken in
$N=2$ and $N=4$ rigid supersymmetric models and in the Wess-Zumino model and
super QED [4,5] for low numbers of loops.

Supersymmetry anomalies do exist and they were first found [6] in the context
of two dimensional models which also broke translation symmetry, that is,
gravitational anomalies [7]. The supersymmetry and gravitational anomalies
occurred in the same anomaly multiplet. A supersymmetry anomaly was found in
four dimensions in reference [8]. However, this appeared in theories that had
a gauge anomaly and as the calculation was performed in the Wess-Zumino gauge,
which is not a supersymmetric gauge choice, the anomaly just reflected the
presence of the gauge anomaly rather than a genuine supersymmetry anomaly.

Recently it has been proposed that when one couples rigid supersymmetric
models to supergravity then one can find a supersymmetry anomaly and this
occurs even for the free Wess-Zumino model coupled to supergravity [9-12].
Unlike the previously mentioned papers this work included supergravity
couplings. Indeed it was claimed that there was an anomaly in the free Wess-
Zumino model in the presence of background supergravity fields [11]. A
supersymmetric theory possesses spacetime translations and, as the presence of
general relativity requires local translations, it follows from the
supersymmetry algebra that supersymmetry must also be a local symmetry. To
have an anomaly in a local symmetry for which there are gauge particles, in
this case the gravitini, renders the theory to be inconsistent. As a result if
the claim of references [9-12] were to be true then the research on
supersymmetry of the past fifty years would have been for nothing. It is
therefore of some interest to see if the claims are true. Very recently
certain one loop current correlation functions of the currents in the Wess-
Zumino were calculated using a Pauli-Villars regulator. It was found that
there was no anomaly in the supersymmetry current [13] in these correlation
functions. The results of references [9-12] were also further discussed in
reference [14].

In this paper we will consider the free Wess-Zumino model in the presence of a
background supergravity background, in fact old minimal supergravity [15]. The
corresponding path integral is

$e^{iZ}=\int{\cal D}\varphi e^{iS^{WZ}+j\cdot h}$ $None$

The action for the free Wess-Zumino model is given by

$S^{WZ}=\int d^{4}x\bigg{\\{}-{1\over
2}\left(\partial_{\mu}A\right)^{2}-{1\over
2}\left(\partial_{\mu}B\right)^{2}-{1\over
2}\bar{\chi}\partial\\!\\!\\!/\chi\bigg{\\}}$ $None$

and the integration ${\cal D}\varphi$ is over the fields in the Wess-Zumino
model, namely $A,B$ and $\chi_{\alpha}$. We do not need the auxiliary fields
$F$ and $G$. We are using the conventions of reference [16] where a full
explanation of many of the equations in this paper can also be found. The term
$j\cdot h$ represents the coupling of the fields of the Wess-Zumino model to
the background supergravity fields and it is given by

$j\cdot h\equiv i\int
d^{4}x(2h^{\mu\nu}\theta_{\mu\nu}+\bar{\psi}^{\mu\alpha}j_{\mu\alpha}-b^{\mu}j^{(5)}_{\mu})$
$None$

Here the supergravity fields $h^{\mu\nu}$, $\psi^{\mu\alpha}$ and $b^{\mu}$
multiply the energy momentum tensor, the supersymmetry current and the axial
current respectively. We do not need the two auxiliary fields $M$ and $N$. For
the free Wess-Zumino model the currents are given by

$\theta_{\mu\nu}={2\over
3}\left(\partial_{\mu}A\partial_{\nu}A+\partial_{\mu}B\partial_{\nu}B\right)+{1\over
4}\bar{\chi}\left(\gamma_{\mu}\partial_{\nu}+\gamma_{\nu}\partial_{\mu}\right)\chi$
$-{1\over
6}\eta_{\mu\nu}\left(\partial_{\rho}A\partial^{\rho}A+\partial_{\rho}B\partial^{\rho}B\right)-{1\over
3}A\partial_{\mu}\partial_{\nu}A-{1\over
3}B\partial_{\mu}\partial_{\nu}B+{1\over
3}\eta_{\mu\nu}\left(A\partial^{2}A+B\partial^{2}B\right)$ $j_{\mu}={1\over
3}\partial\\!\\!\\!/\left(A-i\gamma_{5}B\right)\gamma_{\mu}\chi+{2\over
3}\partial_{\mu}\left(A+i\gamma_{5}B\right)\chi-{2\over
3}\left(A+i\gamma_{5}B\right)\partial_{\mu}\chi+{2\over
3}\left(A+i\gamma_{5}B\right)\gamma_{\mu}\partial\\!\\!\\!/\chi$
$j_{\mu}^{(5)}={2\over 3}\left(B\partial_{\mu}A-A\partial_{\mu}B+{i\over
4}\bar{\chi}\gamma_{5}\gamma_{\mu}\chi\right)$ $None$

It is straight forward to verify that they are conserved in the free Wess-
Zumino model. These currents are known to belong to a supermultiplet [17] and
their supersymmetry transformations are given by

$\delta
j_{\mu}^{(5)}=i\bar{\epsilon}\gamma_{5}j_{\mu},\quad\delta\theta_{\mu\nu}=-{1\over
4}\bar{\epsilon}\left(\gamma_{\mu\rho}\partial^{\rho}j_{\nu}+\gamma_{\nu\rho}\partial^{\rho}j_{\mu}\right),$
$\delta
j_{\mu}=\left(2\gamma^{\nu}\theta_{\mu\nu}-i\gamma_{5}\partial\\!\\!\\!/j_{\mu}^{(5)}-{i\over
2}\gamma_{5}\gamma_{\mu\kappa\lambda}\partial^{\kappa}j^{\lambda(5)}\right)\epsilon$
$None$

The transformations of the old minimal supergravity fields at the linearised
level can be found in that book [16].

At the classical level, the action and background coupling shown in equation
(1) are supersymmetric in four dimensions. However, while this is true outside
four dimensions for the free Wess-Zumino action it is not true for the field
background supergravity coupling. Indeed the supersymmetry transformations of
equation (4) are not true outside four dimensions as they require a Fierz
reshuffle. One instead finds that under a supersymmetry variation

$\delta(j\cdot h)\equiv G=i\int
d^{4}xP^{\mu\alpha\beta}\partial_{\mu}\chi_{\alpha}\chi_{\beta}$ $None$

plus a term that involves $\gamma^{\mu}j_{\mu}$ and one that involves
$\partial^{\mu}j_{\mu}^{(5)}$. We find that

$P^{\mu\alpha\beta}=(-{1\over
3}\bar{\epsilon}^{\alpha}\left(\bar{\psi}^{\lambda}\gamma^{\mu}\gamma_{\lambda}\right)^{\beta}-{2\over
3}\bar{\epsilon}^{\beta}\bar{\psi}^{\mu\alpha}-{2\over
3}\bar{\epsilon}^{\alpha}\bar{\psi}^{\mu\beta})+\big{\\{}\bar{\epsilon}\rightarrow\bar{\epsilon}i\gamma_{5},\bar{\psi}^{\mu}\rightarrow\bar{\psi}^{\mu}i\gamma_{5}\big{\\}}$
$-{1\over
3}\bar{\epsilon}\gamma^{\mu}\psi^{\nu}\left(C\gamma_{\nu}\right)^{\alpha\beta}-{2\over
3}\bar{\epsilon}\gamma^{\nu}\psi^{\mu}\left(C\gamma_{\nu}\right)^{\alpha\beta}-{1\over
3}\bar{\epsilon}\gamma^{\mu}\gamma_{5}\psi^{\nu}\left(C\gamma_{\nu}\gamma_{5}\right)^{\alpha\beta}$
$None$

We now recall [4] the general procedure for calculating in a theory with a
symmetry violating regulator which we take to be dimensional regularisation
but the procedure can be readily applied when other regulators are used. In
particular we take the spacetime dimension $D$ to be $D=4-2\epsilon$. We
consider a theory with fields $\varphi$ and an action $S$ whose variation is
given by $\delta S=G$. The quantity $G$ is only non-zero due to the presence
of the regulator and so it will vanish as the regulator is removed. In our
case as $\epsilon\to 0$. The one particle irreducible effective action
$\Gamma$ derived from the path integral with this action satisfies the Ward
identity

$\delta\varphi{\delta\Gamma\over\delta\varphi}=G\cdot\Gamma$ $None$

where the term on the right corresponds to the insertion of the operator $G$.
To renormalise the theory we consider the quantities

$\Gamma_{ren}=\Gamma-\Gamma_{div}-\Gamma^{FLC},\quad\quad(G\cdot\Gamma)_{ren}=G\cdot\Gamma-(G\cdot\Gamma)_{div}$
$None$

where $\Gamma_{div}$ and $(G\cdot\Gamma)_{div}$ are the divergent parts of
$\Gamma$ and $(G\cdot\Gamma)$ respectively and so they have a one over
$\epsilon$ pole. The term $\Gamma^{FLC}$ allows for the possibility to add to
the effective action a term which is finite and local in the fields. Clearly
the coefficient of the pole does not depend on $\epsilon$ and so we can take
this to be in four dimensions and as such it obeys the Ward identity. As $G$
vanishes as $\epsilon$ goes to zero the term $(G\cdot\Gamma)_{ren}$ vanishes
in this limit. However, $(G\cdot\Gamma)_{div}$ is of the form
${\epsilon\over\epsilon}$ and so this can be finite and non-zero. We will find
that in the limit as $\epsilon$ goes to zero the renormalised effective action
obeys the Ward identity

$\delta\varphi{\delta\Gamma_{ren}\over\delta\varphi}=0$ $None$

provided there do exist finite local counter terms in $\Gamma^{FLC}$ which can
be used to remove any parts of $(G\cdot\Gamma)_{div}$ which are not invariant
under the symmetry. If this is not the case then the theory has an anomaly.

We now apply this general theory to supersymmetry and the free Wess-Zumino
model in the background of supergravity fields reviewed above. In particular
we will consider the terms in $(G\cdot\Gamma)_{div}$ which depend on the
auxiliary field $b_{\mu}$. They also depend on the supersymmetry parameter
$\epsilon^{\alpha}$ and the gravitino $\psi^{\mu}_{\alpha}$ as these appear in
the expression for $G$ of equation (7). The corresponding lowest order such
Feynman diagram has two vertices, one contains the background gravitino and
the two spinors $\chi_{\alpha}$ which form internal propagators. The other
vertex contains the coupling $b^{\mu}j^{(5)}_{\mu}$ that appears in the
original action of equation (1). The chiral current $j^{(5)}_{\mu}$ given in
equation (4) contains the same two chiral spinors which also form propagators.
The resulting Feynman graph corresponds to the expression

${1\over 3}\int d^{4}p\int{d^{4}k\over(2\pi)^{4}}{1\over
k^{2}(p+k)^{2}}Tr\\{{{k}}\\!\\!\\!/k_{\mu}CP^{\mu}(p)(p\\!\\!\\!/+k\\!\\!\\!/)\gamma_{5}\gamma_{\nu}\\}b^{\nu}(-p)$
$None$

where $C$ is the charge conjugation matrix.

The graph of equation (11) is relatively straight forward to evaluate using
for example the formula in the classic Diagrammer [18] or the book [19]. One
finds that it is equal to

$-{2\over 9}\int d^{4}p{1\over(4\pi)^{2}}{p^{2}\over
12}\bigg{\\{}-3\bar{\epsilon}\gamma_{5}b\\!\\!\\!/p\\!\\!\\!/\gamma\cdot\psi+2\bar{\epsilon}\gamma_{5}\gamma\cdot\psi
b\cdot
p+2\bar{\epsilon}\gamma_{5}b\\!\\!\\!/p\cdot\psi-2\bar{\epsilon}\gamma_{5}p\\!\\!\\!/b\cdot\psi\bigg{\\}}$
$None$

This can be removed by the finite local counter term

$\Gamma^{FLC}={2\over 9}\int d^{4}p{1\over(4\pi)^{2}}{p^{2}\over
12}i\bigg{\\{}-\bar{\psi}_{\mu}p\\!\\!\\!/\psi^{\mu}+2\bar{\psi}\cdot
p\gamma\cdot\psi+{5\over 2}\bar{\psi}\cdot\gamma
p\\!\\!\\!/\gamma\cdot\psi\bigg{\\}}$ $None$

The two additional terms mentioned below equation (3) do not lead to any
contribution to the effective action and this corresponds to the fact that
S-supersymmetry and the chiral current are conserved at least for the
processes being considered in this paper. As such there is no supersymmetry
anomaly in this contribution to the effective action.

The calculation in this paper is only a first step in showing the consistency
of supersymmetry in a supergravity background using dimensional
regularisation. It would be interesting to calculate the contribution
containing the graviton. Indeed as there are no gauge fields which require
gauge fixing this calculation and the one in this paper could be
simultaneously preformed in superspace using super Feynman rules. It may also
be possible to completely evaluate the path integral of equation (1) as the
contributions are only at one loop.

If one has gauge particles and one wants to calculate in components of
familiar form then one must fix the Wess-Zumino gauge. However, this
explicitly breaks supersymmetry. How to take account of this in the Ward
Identities was explained in references [4] and [5] and was successfully
implemented in the context of the Wess-Zumino model, super QED and extended
supersymmetry gauge theories although this was in the absence of background
supergravity. It would be interesting to compare this approach to this problem
with that of reference [20].

In the paper [6] it was stated concerning in four dimensions that ”if the
anomalies lie in a supermultiplet and the multiplet of currents is contained
in a vector superfield then there are no supersymmetry anomalies. The proof of
this theorem examines the most general vector multiplet that has a chiral
current as its first component. One then identifies the most general possible
energy-momentum tensor and demands that it be symmetric and conserved. The
supervariation of these constraints imply that there exists a conserved
supercurrent. ” This statement would appear to transcend the theory being
considered and it would mean that there are no supersymmetry anomalies in four
dimensions as there are no gravitational anomalies. It would be good to
consider this in detail in the presence of background supergravity. A review
of possible anomaly multiplets can be found in chapter twenty of reference
[16] although this was before the discovery of supersymmetry anomalies.

Acknowledgments I wish to thank the SFTC for support from Consolidated grants
number ST/T000759/1 and ST/P000258/1.

References

[1] J. Ashmore, A Method of Gauge Invariant Regularization, Nuovo Cimento
Lett. 4 (1972) 289.

[2] C. Bollini and J. Giambiagi, Dimensional Renormalization: The Number of
Dimensions as a Regularizing Parameter, Nuovo Cimento 12B (1972) 20;

[3] I. Jack, Two Loop Beta Functions for Supersymmetric Gauge Theories, Phys.
Lett. 147B (1984) 405

[4] P.S. Howe, A.J. Parkes and P.C. West, Dimensional Regularisation and
Supersymmetry, Phys. Lett. 147B (1984) 409.

[5] P.S. Howe, A.J. Parkes and P.C. West, On Supersymmetry Anomalies, Phys.
Lett. 150B (1985) 149.

[6] P. Howe and P. West, Gravitational Anomalies in Supersymmetric Theories,
Phys. Lett. 156B 335 (1985).

[7] L. Alvarez-Gaume and E. Witten, Gravitational anomalies, Nucl. Phy B234
(1984) 269.

[8] H. Itoyama, V. P. Nair, and H.-c. Ren, Nucl. Phys.Supersymmetry Anomalies
and some aspects of renormalisation, B262, 317 (1985).

[9] I. Papadimitriou, Supercurrent anomalies in 4d SCFTs, JHEP07(2017) 038,
ArXiv:1703.04299.

[10] I. Papadimitriou, Supersymmetry anomalies in new minimal supergravity,
JHEP09(2019)039 [1904.00347].

[11] G. Katsianis, I. Papadimitriou, K. Skenderis and M. Taylor, Anomalous
Supersymmetry, Phys. Rev. Lett.122(2019) 231602, ArXiv:1902.06715.

[12] I. Papadimitriou, Supersymmetry anomalies inN= 1conformal supergravity,
JHEP04(2019) 040. ArXiv:1902.06717.

[13] A. Bzowski, G. Festuccia and V. Prochazka. Consistency of supersymmetric
t Hooft anomalies, ArXiv:2011.09978

[14] G. Katsianis, I. Papadimitriou, K. Skenderis and M. Taylor, Supersymmetry
anomaly in the superconformalWess-Zumino model

[15] K.S. Stelle and P. West. Minimal Auxiliary Fields for Supergravity, Phys.
Lett 74B, 330 (1978); S. Ferrara and P. van Nieuwenhuizen, The Auxiliary
Fields of Supergravity, Phys. Lett. B74, 333 (1978).

[16] P. West, Introduction to Supersymmetry and Supergravity, World Scientific
Publishing 1986, second edition 1990.

[17] S. Ferrara and B. Zumino, Transformation Properties of the Supercurrent,
Nucl. Phys. B87(1975) 207.

[18] G. ’t Hooft and M. Veltman, Diagrammer, CERN preprint 73-9, 3 September
1973.

[19] P. Ramond, Field Theory: A Modern Primer, Benjamin/Cummings 1981.

[20] S.M. Kuzenko, A. Schwimmer and S. Theisen, Comments on Anomalies in
SupersymmetricTheories, J. Phys. A53 (2020) 064003, ArXiv:1909.07084.


# Reduction of back switching by large damping ferromagnetic material

Tomohiro Taniguchi1, Yohei Shiokawa2, and Tomoyuki Sasaki2 1National Institute
of Advanced Industrial Science and Technology (AIST), Research Center for
Emerging Computing Technologies, Tsukuba, Ibaraki 305-8568, Japan,
2TDK Corporation, Advanced Products Development Center, Ichikawa, Chiba
272-8558, Japan

###### Abstract

Recent studies on magnetization dynamics induced by spin-orbit torque have
revealed a weak dependence of the critical current for magnetization switching
on the damping constant of a ferromagnetic free layer. This study, however,
reveals that the damping constant nevertheless plays a key role in
magnetization switching induced by spin-orbit torque. An undesirable
switching, returning to an initial state, named as back switching, occurs in a
ferromagnet with an easy axis parallel to the current direction. Numerical and
theoretical analyses reveal that back switching is strongly suppressed when
the damping constant of the ferromagnet is large.

Low-damping ferromagnetic materials has been investigated for spintronics
applications Oogane et al. (2006); Bilzer et al. (2006); Mizukami et al.
(2008); Iihama et al. (2012); Konoto et al. (2013); Tsunegi et al. (2014).
These materials are interesting because the critical current for exciting
magnetization dynamics by spin-transfer torque Slonczewski (1996); Berger
(1996) in two-terminal devices is typically proportional to the damping
constant Sun (2000); Grollier et al. (2003), and therefore, low-damping
ferromagnetic materials help to reduce power consumption. The value of the
damping constant in typical ferromagnets, such as Fe, Co, Ni, and their
alloys, is on the order of $10^{-3}-10^{-2}$, even after the effect of spin
pumping Tserkovnyak et al. (2002); Mizukami et al. (2002) is taken into
account Konoto et al. (2013); Tsunegi et al. (2014). Studies on three-terminal
devices manipulated by spin-orbit torque Liu et al. (2012); Pai et al. (2012);
Yu et al. (2014); Cubukcu et al. (2014); You et al. (2015); Torrejon et al.
(2015); Fukami et al. (2016a, b); Lau et al. (2016); van den Brink et al.
(2016); Oh et al. (2016); Shiokawa et al. (2019), however, have revealed that
the dependence of the critical current on the damping constant is weak when
the easy axis of the ferromagnet is perpendicular to the film plane Lee et al.
(2013); Taniguchi (2019) or parallel to the current direction Taniguchi
(2020); these systems are named type Z and type X, respectively, in Ref.
Fukami et al. (2016a). Then, a question arises as to what the role of the
damping constant is in three-terminal devices.

The purpose of this study is to investigate the relation between the
magnetization state manipulated by the spin-orbit torque and the magnetic
damping constant by solving the Landau-Lifshitz-Gilbert (LLG) equation. We
focus on the type-X system because the switching mechanism in the system is
not fully understood yet, unlike the other systems Sun (2000); Grollier et al.
(2003); Lee et al. (2013); Taniguchi (2019). An undesirable return of the
magnetization, called back switching in this paper, occurs in a ferromagnet,
where, although the spin-orbit torque brings the magnetization close to the
switched state, it returns to the initial state after turning off the current.
As a result, the phase diagram of the magnetization state as a function of the
current and external magnetic field alternately shows switched and non-
switched states. It is revealed that the back switching appears as a result of
the magnetization precession around the perpendicular axis after the current
is turned off. The back-switching region is strongly suppressed by using
large-damping ferromagnetic materials.

Figure 1:  (a) Schematic illustration of the system. The spin Hall effect in
the bottom nonmagnet generates spin current with a polarization along the $y$
direction and excites the magnetization ($\mathbf{m}$) dynamics in the top
ferromagnet. The easy axis of the ferromagnet is parallel to the $x$ axis. An
external magnetic field $H_{\rm appl}$ is applied in the $z$ direction. (b)
Example of the time evolution of $m_{x}$ for $j=300$ MA cm-2 and $H_{\rm
appl}=500$ Oe.

The system studied in the work is an in-plane magnetized ferromagnet placed on
a nonmagnetic metal and schematically shown in Fig. 1(a). Electric current
flowing in the nonmagnet in the $x$ direction generates pure spin current that
is injected into the ferromagnet with the spin polarization in the $y$
direction. The spin current excites the spin-transfer torque acting on the
magnetization in the ferromagnet and induces magnetization dynamics. The $z$
axis is perpendicular to the film plane. We denote the unit vectors pointing
in the magnetization direction in the ferromagnet and along the $k$-axis
($k=x,y,z$) as $\mathbf{m}$ and $\mathbf{e}_{k}$, respectively. The
magnetization dynamics are described by the LLG equation

$\frac{d\mathbf{m}}{dt}=-\gamma\mathbf{m}\times\mathbf{H}-\gamma H_{\rm
s}\mathbf{m}\times\left(\mathbf{e}_{y}\times\mathbf{m}\right)+\alpha\mathbf{m}\times\frac{d\mathbf{m}}{dt},$
(1)

where $\gamma$ and $\alpha$ are the gyromagnetic ratio and the Gilbert damping
constant, respectively. The magnetic field $\mathbf{H}$ consists of the in-
plane magnetic anisotropy field $H_{\rm K}$ in the $x$ direction, the
demagnetization field $-4\pi M$ in the $z$ direction, and an external field
$H_{\rm appl}$ in the $z$ direction:

$\mathbf{H}=H_{\rm K}m_{x}\mathbf{e}_{x}+\left(H_{\rm appl}-4\pi
Mm_{z}\right)\mathbf{e}_{z}.$ (2)

The strength of the spin-transfer torque is given by

$H_{\rm s}=\frac{\hbar\vartheta j}{2eMd},$ (3)

where $\vartheta$ is the spin Hall angle in the nonmagnet, whereas $M$ and $d$
are the saturation magnetization and thickness of the ferromagnet. The
electric current density is denoted as $j$. The values of the parameters used
in this study were taken from typical experiments Liu et al. (2012); Pai et
al. (2012); Yu et al. (2014); Cubukcu et al. (2014); You et al. (2015);
Torrejon et al. (2015); Fukami et al. (2016a, b); Lau et al. (2016); van den
Brink et al. (2016); Oh et al. (2016) as $M=1500$ emu c.c.-1, $H_{\rm K}=200$
Oe, $\gamma=1.764\times 10^{7}$ rad Oe-1 s-1, $\vartheta=0.4$, and $d=1.0$ nm.

The magnetic field $\mathbf{H}$ is related to the magnetic energy density $E$
via $E=-M\int d\mathbf{m}\cdot\mathbf{H}$,

$E=-MH_{\rm appl}m_{z}-\frac{MH_{\rm K}}{2}m_{x}^{2}+2\pi M^{2}m_{z}^{2}.$ (4)

The energy density $E$ has two minima at
$\mathbf{m}_{0\pm}=\pm\sqrt{1-m_{0z}^{2}}\mathbf{e}_{x}+m_{0z}\mathbf{e}_{z}$,
where $m_{0z}=H_{\rm appl}/(H_{\rm K}+4\pi M)$. Throughout this paper, the
initial state is set to be $\mathbf{m}_{0-}$, which points in the negative $x$
direction. Accordingly, we call the other stable state, $\mathbf{m}_{0+}$
pointing in the positive $x$ direction, the switched state. Thus, we are
interested in experiments where the initial state is reset to
$\mathbf{m}_{0-}$ during each trial of magnetization switching. By convention,
we will focus on switching by a positive current and field region.

Figure 2:  Phase diagrams of $m_{x}$ (a) at a fixed point in the presence of
current and (b) in a relaxed state after turning off the current. The initial
state is $\mathbf{m}_{0-}$. The white dotted line in (a) represents the
theoretical formula of the critical current density derived in Ref. Taniguchi
(2020). Switching regions are distinguished by the labels Sn ($n=0,1,2$) in
(b). The damping constant $\alpha$ is $0.005$.

Figure 1(b) shows an example of the time evolution of $m_{x}$ in the presence
of current. The damping constant is $\alpha=0.005$. The results indicate that
the magnetization saturates to a fixed point. Figure 2(a) is a phase diagram
summarizing the fixed point of $m_{x}$ in the presence of current, where the
vertical and horizontal axes represent the current density $j$ and the
external magnetic field $H_{\rm appl}$. The magnetization stays near the
initial state [$m_{x}\simeq-1$ shown in black in Fig. 2(a)] in a relatively
small current region, where its boundary is well explained by the critical-
current formula $j_{\rm c}$ Taniguchi (2020) shown by the white dotted line in
Fig. 2(a). On the other hand, the magnetization state above the critical
current satisfies $m_{x}>0$. Therefore, one might suppose that the
magnetization relaxes to the switched state, $m_{x}\simeq+1$, after turning
off the current.

However, the phase diagram after turning off the current shown in Fig. 2(b)
reveals a complicated dependence of the relaxed state on the current and
external magnetic field. There are stripes distinguishing the switched
($\mathbf{m}_{0+}$, yellow) and non-switched ($\mathbf{m}_{0-}$, black)
states. This result indicates that the magnetization returns to the initial
state under certain conditions. We name this phenomenon back switching and
investigate its relation to the damping constant in the following.

Let us briefly comment on another phenomenon, called backhopping Sun et al.
(2009); Oh et al. (2009); Skowroński et al. (2013); Abert et al. (2018);
Safranski and Sun (2019), to avoid any confusion. Backhopping is a phenomenon
in which, after magnetization switching in a free layer is achieved, an
undesirable return to the original resistance-state occurs at a high-bias
voltage. There are two differences between back switching and backhopping.
First, backhopping in two-terminal devices originates from magnetization
switching in the reference layer. On the other hand, back switching relates to
the magnetization dynamics in the free layer. Second, backhopping has often
been investigated by sweeping the current. On the other hand, we reset the
initial state of the magnetization at each trial because we are interested in
designing the operation conditions of memory devices. Back switching
originates from the fact that the fixed point of the magnetization in the
presence of current is not parallel to the easy axis of the magnetization, as
discussed below, and thus, it will occur in not only type-X devices, but also
type-Z devices Lee et al. (2013); Taniguchi (2019).

Note that the horizontal stripe in the relatively low current region of Fig.
2(b), slightly above $j_{\rm c}$, was analyzed in Ref. Taniguchi (2020), and
therefore, it will be excluded from the following discussion. We label the
other switched regions shown by the vertical stripes as the Sn-region
($n=0,1,2,\cdots$); see Fig. 2(b). The role of the integer $n$ will be
clarified below.

Figure 3:  (a) Dynamic trajectory in the S0-region; $H_{\rm appl}=500$ Oe, and
$j=300$ MA cm-2. The red line represents the trajectory from the initial state
to the fixed point (yellow circle) in the presence of the current, whereas the
blue line corresponds to the relaxation dynamics after turning off the
current. The symbols $\mathbf{m}_{0\pm}$ indicate the locations of the initial
($\mathbf{m}_{0-}$) and switched ($\mathbf{m}_{0+}$) states. (b) Time
evolutions of $m_{x}$ (black solid line) and normalized energy density
$\varepsilon$ (dashed red line) after turning off the current. The red
triangle indicates the time at which $\varepsilon$ becomes zero. (c) Dynamic
trajectory and (d) time evolution of $m_{x}$ and $\varepsilon$ in the region
sandwiched by the S0\- and S1-regions; $H_{\rm appl}=750$ Oe, and $j=300$ MA
cm-2. (e) Dynamic trajectory and (f) time evolution of $m_{x}$ and
$\varepsilon$ in the S1-region; $H_{\rm appl}=1000$ Oe, and $j=300$ MA cm-2.

Figure 3(a) shows the dynamic trajectory of the magnetization obtained by
numerically solving Eq. (1), where $j=300$ MA cm-2 and $H_{\rm appl}=500$ Oe
correspond to the S0 region. The red line represents the trajectory from the
initial state to the fixed point in the presence of the current, whereas the
blue line shows the relaxation dynamics after turning off the current. The
fixed point satisfying $d\mathbf{m}/dt=\bm{0}$ in the presence of the current
is indicated by a yellow circle. Note that the fixed point is far away from
the energetically stable states, $\mathbf{m}_{0\pm}$. This is because the
spin-transfer torque drives the magnetization in the $y$ direction. The black
solid line in Fig. 3(b) shows the time evolution of $m_{x}$ after turning off
the current. Starting from the fixed point located in the region of $m_{x}>0$,
the magnetization relaxes to the switched state $\mathbf{m}_{0+}$ after
showing several precessions around it. Figures 3(c) and 3(d), on the other
hand, show an example of back switching, where the parameter $H_{\rm
appl}=750$ Oe corresponds to the region sandwiched by the S0 and S1-regions.
In this case, the magnetization after turning off the current shows a
precession around the $z$ axis. As a result, even though the fixed point in
the presence of the current is in the region of $m_{x}>0$, the magnetization
returns to the initial state ($\mathbf{m}_{0-}$). Let us then move to the next
switching region, the S1-region. The dynamic trajectory, as well as $m_{x}$,
shown in Figs. 3(e) and 3(f), indicates that, starting from the fixed point in
the region of $m_{x}>0$, the magnetization returns once to the region of
$m_{x}<0$ before relaxing to the switched state. Now the meaning of the
integer $n$ ($n=0,1,2,\cdots$) we used to distinguish the switched regions
becomes clear: it represents how many times the magnetization returns to the
region of $m_{x}<0$ before relaxing to the switched state.

These results imply that the precession around the $z$ axis after turning off
the current is the origin of back switching. Such a precession is induced by
the precession torque, $-\gamma\mathbf{m}\times\mathbf{H}$. In particular, the
precession around the $z$ axis occurs when the energy density at the fixed
point of the magnetization in the presence of the current is larger than the
saddle-point energy given by

$E_{\rm d}=-\frac{MH_{\rm appl}^{2}}{8\pi M}.$ (5)

In fact, a fixed point of the LLG equation is given by

$\displaystyle m_{x}=\frac{H_{\rm appl}H_{\rm s}}{H_{\rm s}^{2}-H_{\rm K}4\pi
M},$ $\displaystyle m_{z}=-\frac{H_{\rm appl}H_{\rm K}}{H_{\rm s}^{2}-H_{\rm
K}4\pi M}.$ (6)

Substituting Eq. (6) into Eq. (4), the energy density at this fixed point is

$E_{j}=\frac{MH_{\rm appl}^{2}H_{\rm K}}{2(H_{\rm s}^{2}-H_{\rm K}4\pi M)}.$
(7)

Note that $H_{\rm s}^{2}-H_{\rm K}4\pi M>0$ because the back switching appears
in the current region above $j_{\rm c}$ for $H_{\rm appl}=0$, which is given
by $j_{\rm c}=[2eMd/(\hbar\vartheta)]\sqrt{H_{\rm K}4\pi M}$ Taniguchi (2020).
Therefore, Eq. (7) is always positive, whereas the saddle-point energy density
given by Eq. (5) is negative. Accordingly, the fixed point in the presence of
the current is in the unstable region corresponding to an energy larger than
the saddle-point energy. As a result, after turning off the current, the
magnetization starts to precess around the $z$ axis, as mentioned above.
Simultaneously, the magnetization loses energy due to damping torque, where
the dissipated energy density is given by

$\Delta E(t)=\frac{\alpha\gamma
M}{1+\alpha^{2}}\int_{0}^{t}dt^{\prime}\left[\mathbf{H}^{2}-\left(\mathbf{m}\cdot\mathbf{H}\right)^{2}\right].$
(8)

Here, $t$ is the time after the current is turned off. When the condition,

$E_{j}-\Delta E(t)=E_{\rm d},$ (9)

is satisfied, the magnetization relaxes to the nearest stable state. The point
here is that the magnetization alternately comes close to two stable states,
$\mathbf{m}_{0\pm}$, due to the precession around the $z$ axis before Eq. (9)
is satisfied. As a result, both $\mathbf{m}_{0+}$ and $\mathbf{m}_{0-}$ can be
the final state, which leads to back switching. The conclusions of the
discussion are confirmed by the dashed red lines in Figs. 3(b), 3(d), and
3(f), where are the time evolutions of the energy density. Here, we introduce
the normalized energy density,

$\varepsilon(t)=\frac{E(t)-E_{\rm d}}{4\pi M^{2}}.$ (10)

The energy density $E$ is a function of time with the initial condition of
$E(t=0)=E_{j}$ because $E$, given by Eq. (4), depends on the magnetization
direction, and the magnetization changes direction in accordance with the LLG
equation. The condition given by Eq. (9) is $\varepsilon(t)=0$ in terms of the
normalized energy density. The times at which $\varepsilon$ becomes zero are
indicated by the red triangles in the figures. It is clear from Figs. 3(b),
3(d), and 3(f) that the final state of the magnetization is determined by
whether $m_{x}$ is positive or negative when $\varepsilon$ becomes zero.

Figure 4:  Time evolutions of $m_{x}$ (black solid lines) and normalized
energy density $\varepsilon$ (red dashed lines) after turning off the current,
where $j=300$ MA cm-2 and (a) $H_{\rm appl}=750$ Oe and (b) $1000$ Oe. The
damping constant $\alpha$ is $0.050$. The red triangles indicate the time at
which $\varepsilon$ becomes zero. Phase diagrams of $m_{x}$ (c) at a fixed
point in the presence of current and (d) in a relaxed state after turning off
the current.

It is expected that using a large-damping ferromagnetic material will result
in a reduction of the back-switching region. Remember that back switching
occurs due to the precession around the $z$ axis. When the damping constant is
large, the energy dissipation given by Eq. (8) becomes large within a short
time $t$, and the condition given by Eq. (9) is immediately satisfied before
the magnetization returns to the region of $m_{x}<0$ by the precession around
the $z$ axis. Accordingly, back switching does not occur. On the other hand,
when the damping constant is small, it takes a long time to dissipate the
energy to satisfy Eq. (9), during which time the magnetization shows the
precession around the $z$ axis. As a result, back switching appears.

To verify this picture, we evaluated the phase diagram of the magnetization
state for a relatively large damping constant, $\alpha=0.050$, which is ten
times large than that used above. Figures 4(a) and 4(b) show the time
evolutions of $m_{x}$ and $\varepsilon$ after turning off the current, where
$H_{\rm appl}$ is (a) $750$ Oe and (b) $1000$ Oe. Contrary to the dynamics
shown in Figs. 3(d) and 3(f), the magnetization in Figs. 4(a) and 4(b) does
not return to the region of $m_{x}<0$ because the energy dissipates quickly
due to the large damping torque. As a result, the magnetization immediately
relaxes to the switched state. Figures 4(c) and 4(d) summarize the
magnetization state in the presence of current and after turning off the
current, respectively. Comparing Fig. 4(c) with Fig. 2(a), it is clear that
the phase diagram of the magnetization state in the presence of the current is
approximately independent of the damping constant. On the other hand, the
back-switching region is strongly suppressed, as can be seen by the comparing
Fig. 2(b) and 4(d). This is because the large damping torque immediately
dissipates the energy from the ferromagnet and leads to a fast relaxation to
the switched state, due to which the S0-region dominates in the phase diagram.

The existence of the back switching gives an upper limit of a write-current
margin for memory applications, and therefore, it restricts the device design
and/or manipulation conditions. On the other hand, back switching might be
applicable to other devices such as a random number generator. The analysis in
this study provides fruitful insights for the development of spintronics
applications utilizing spin-orbit torque.

In conclusion, the phase diagram of the magnetization state in a type-X spin-
orbit torque device was calculated as a function of the electric current
density and the external magnetic field. The magnetization state after turning
off the current has stripe structures alternately showing switched and non-
switched states. Such a non-switched state, named back switching here, occurs
as a result of the magnetization precession around the perpendicular axis
after the current is turned off. The back-switching region is reduced by using
large-damping ferromagnetic materials, whereas the critical current inducing
the magnetization switching is approximately independent of the damping
constant.

The authors thank to Masamitsu Hayashi, Shinji Isogami, Toru Oikawa, Takehiko
Yorozu, and Seiji Mitani for valuable discussion. This work was supported by
funding from TDK Corporation.

## References

  * Oogane et al. (2006) M. Oogane, T. Wakitani, S. Yakata, R. Yilgin, Y. Ando, A. Sakuma, and T. Miyazaki, Jpn. J. Appl. Phys. 45, 3889 (2006).
  * Bilzer et al. (2006) C. Bilzer, T. Devolder, J.-V. Kim, G. Counil, C. Chappert, S. Cardoso, and P. P. Freitas, J. Appl. Phys. 100, 053903 (2006).
  * Mizukami et al. (2008) S. Mizukami, H. Abe, D. Watanabe, M. Oogane, Y. Ando, and T. Miyazaki, Appl. Phys. Express 1, 121301 (2008).
  * Iihama et al. (2012) S. Iihama, Q. Ma, T. Kubota, S. Mizukami, Y. Ando, and T. Miyazaki, Appl. Phys. Express 5, 083001 (2012).
  * Konoto et al. (2013) M. Konoto, H. Imamura, T. Taniguchi, K. Yakushiji, H. Kubota, A. Fukushima, K. Ando, and S. Yuasa, Appl. Phys. Express 6, 073002 (2013).
  * Tsunegi et al. (2014) S. Tsunegi, H. Kubota, S. Tamaru, K. Yakushiji, M. Konoto, A. Fukushima, T. Taniguchi, H. Arai, H. Imamura, and S. Yuasa, Appl. Phys. Express 7, 033004 (2014).
  * Slonczewski (1996) J. C. Slonczewski, J. Magn. Magn. Mater. 159, L1 (1996).
  * Berger (1996) L. Berger, Phys. Rev. B 54, 9353 (1996).
  * Sun (2000) J. Z. Sun, Phys. Rev. B 62, 570 (2000).
  * Grollier et al. (2003) J. Grollier, V. Cros, H. Jaffrés, A. Hamzic, J. M. George, G. Faini, J. B. Youssef, H. LeGall, and A. Fert, Phys. Rev. B 67, 174402 (2003).
  * Tserkovnyak et al. (2002) Y. Tserkovnyak, A. Brataas, and G. E. W. Bauer, Phys. Rev. Lett. 88, 117601 (2002).
  * Mizukami et al. (2002) S. Mizukami, Y. Ando, and T. Miyazaki, J. Magn. Magn. Mater. 239, 42 (2002).
  * Liu et al. (2012) L. Liu, O. J. Lee, T. J. Gudmundsen, D. C. Ralph, and R. A. Buhrman, Phys. Rev. Lett. 109, 096602 (2012).
  * Pai et al. (2012) C.-F. Pai, L. Liu, Y. Li, H. W. Tseng, D. C. Ralph, and R. A. Buhrman, Appl. Phys. Lett. 101, 122404 (2012).
  * Yu et al. (2014) G. Yu, P. Upadhyaya, Y. Fan, J. G. Alzate, W. Jiang, K. L. Wong, S. Takei, S. A. Bender, L.-T. Chang, Y. Jiang, et al., Nat. Nanotechnol. 9, 548 (2014).
  * Cubukcu et al. (2014) M. Cubukcu, O. Boulle, M. Drouard, K. Garello, C. O. Avci, I. M. Miron, J. Langer, B. Ocker, P. Gambardella, and G. Gaudin, Appl. Phys. Lett. 104, 042406 (2014).
  * You et al. (2015) L. You, O. L. nd D. Bhowmik, D. Labanowski, J. Hong, J. Bokor, and S. Salahuddin, Proc. Natl. Acad. Sci. 112, 10310 (2015).
  * Torrejon et al. (2015) J. Torrejon, F. Garcia-Sanchez, T. Taniguchi, J. Shinha, S. Mitani, J.-V. Kim, and M. Hayashi, Phys. Rev. B 91, 214434 (2015).
  * Fukami et al. (2016a) S. Fukami, T. Anekawa, C. Zhang, and H. Ohno, Nat. Nanotechnol. 11, 621 (2016a).
  * Fukami et al. (2016b) S. Fukami, T. Anekawa, C. Zhang, and H. Ohno, Nat. Mater. 15, 535 (2016b).
  * Lau et al. (2016) Y.-C. Lau, D. Betto, K. Rode, J. M. D. Coey, and P. Stamenov, Nat. Nanotechnol. 11, 758 (2016).
  * van den Brink et al. (2016) A. van den Brink, G. Vermijs, A. Solignac, J. Koo, J. T. Kohlhepp, H. J. M. Swagten, and B. Koopmans, Nat. Commun. 7, 10854 (2016).
  * Oh et al. (2016) Y.-W. Oh, S. h. C. Baek, Y. M. Kim, H. Y. Lee, K.-D. Lee, C.-G. Yang, E.-S. Park, K.-S. Lee, K.-W. Kim, G. Go, et al., Nat. Nanotechnol. 11, 878 (2016).
  * Shiokawa et al. (2019) Y. Shiokawa, E. Komura, Y. Ishitani, A. Tsumita, K. Suda, Y. Kakinuma, and T. Sasaki, AIP Adv. 9, 035236 (2019).
  * Lee et al. (2013) K.-S. Lee, S.-W. Lee, B.-C. Min, and K.-J. Lee, Appl. Phys. Lett. 102, 112410 (2013).
  * Taniguchi (2019) T. Taniguchi, Phys. Rev. B 100, 174419 (2019).
  * Taniguchi (2020) T. Taniguchi, Phys. Rev. B 102, 104435 (2020).
  * Sun et al. (2009) J. Z. Sun, M. C. Gaidis, G. Hu, E. J. O’Sullivan, S. L. Brown, J. J. Nowak, P. L. Trouilloud, and D. C. Worledge, J. Appl. Phys. 105, 07D109 (2009).
  * Oh et al. (2009) S.-C. Oh, S.-Y. Park, A. Manchon, M. Chshiev, J.-H. Han, H.-W. Lee, J.-E. Lee, K.-T. Nam, Y. Jo, Y.-C. Kong, et al., Nat. Phys. 5, 898 (2009).
  * Skowroński et al. (2013) W. Skowroński, P. Ogrodnik, J. Wrona, T. Stobiecki, R. Świrkowicz, J. Barnaś, G. Reiss, and S. van Dijken, J. Appl. Phys. 114, 233905 (2013).
  * Abert et al. (2018) C. Abert, H. Sepehri-Amin, F. Bruckner, C. Vogler, M. Hayashi, and D. Suess, Phys. Rev. Applied 9, 054010 (2018).
  * Safranski and Sun (2019) C. Safranski and J. Z. Sun, Phys. Rev. B 100, 014435 (2019).


# Precise nanoscale temperature mapping in operational microelectronic devices
by  use of a phase change material

Qilong Cheng1111Equal contribution. Sukumar Rajauria21 email:
<EMAIL_ADDRESS>Erhard Schreck2 Robert Smith2 Na Wang2 Jim Reiner2
Qing Dai2 David Bogy1 1Department of Mechanical Engineering, UC Berkeley,
California 94720 USA.  2Western Digital Corporation, Recording Sub System
Staging and Research, San Jose, CA 95135 USA.

###### Abstract

The microelectronics industry is pushing the fundamental limit on the physical
size of individual  elements to produce faster and more powerful integrated
chips. These chips have nanoscale features that dissipate power resulting in
nanoscale hotspots  leading to device failures. To understand the reliability
impact of the hotspots, the device needs to be tested under the actual
operating conditions.  Therefore, the development of high-resolution
thermometry techniques is required to understand the heat dissipation
processes during the device operation. Recently, several thermometry
techniques have been proposed, such as radiation thermometry, thermocouple
based contact thermometry, scanning thermal microscopy (SThM), scanning
transmission electron microscopy (STEM) and transition based threshold
thermometers. However, most of these techniques have limitations including the
need for extensive calibration, perturbation of the actual device temperature,
low throughput, and the use of ultra-high vacuum. Here, we present a facile
technique, which uses a thin film contact thermometer based on the phase
change material $Ge_{2}Sb_{2}Te_{5}$,  to precisely map thermal contours from
the nanoscale to the microscale. $Ge_{2}Sb_{2}Te_{5}$ undergoes a crystalline
transition at Tg with  large changes in its electric conductivity, optical
reflectivity and density. Using this approach, we map the  surface temperature
of a nanowire and an embedded micro-heater on the same chip where  the scales
of the temperature contours differ by three orders of magnitude. The spatial
resolution can be as high as 20 nanometers thanks to the continuous nature of
the thin film.

Nanoscale hotspots, Thermometry, Micro-electronics, Temperature map

The fundamental understanding of thermal dissipation in an integrated chip
HaenschIBM06 ; CavinJNR06 ; PopNR10 ; DejeneNP13 ; JinNL11  requires the
development of a versatile technique capable of reliably mapping the  areal
temperature of various components integrated in the chip ranging from
nanometer to micrometer dimensions BritesNanoscale12 ; LeeSR07 . Various
thermometers were developed to achieve this goal WildeN06 ; DeshpandePRL09 ;
TessierAPL07 ; BritesNanoscale12 ; NonnenmacherAPL92 ; ShiAPL92 ;
MecklenburgScience15 ; BrintlingerNL08 ; GuoNC14 ; GaoN02 ; MartinekU15 ;
MengesNC16 ; KinkhabwalaIEEE15 ; KimACS12 ; HwangRSI14  and can be broadly
classified into two categories: non-contact and contact  based thermometers.
Techniques such as Raman DeshpandePRL09 , fluorescence OkabeNC12 ,
luminescence KucskoNature13  and transmission electron microscopy
MecklenburgScience15 ; IdroboPRL18 ; FengNL18  are non-contact thermometers.
However, the areal resolutions of these methods are limited either by the
optical diffraction limit ChristoffersonJEP08  or by the use of specific
metals and semiconductors MecklenburgScience15 . A scanning thermal microscope
is an extensively used contact thermometer, but it typically suffers from
contact-related artifacts that lead to an under prediction of the device
temperature. This is due to  the thermal coupling strength between the device
and the SThM tip, which is material dependent and difficult to measure
MengesNC16 ; MengesNL12 .

Figure 1: Heat sources inside the head of the hard disk drive. (a) Cross-
sectional cartoon of the head structure showing the embedded heat sources: the
nano-heater and the micro-heater. (b) AFM image of the device: the micro-
heater is embedded and cannot be seen from the surface, while the nano-heater
is located at the center. The dimension of the nano-heater is 1 $\mu$m
$\times$ 20 nm. (c), (d) Simulation: Temperature map of the nano-heater and
the micro-heater with similar peak surface temperature.

Here a novel technique, which uses a phase change material to map the
temperature of an operational microelectronic device, is presented. It
requires minimal effort in temperature calibration and  the temperature
contour  can be mapped using both  contact and non-contact modes such as AFM,
SEM or optical microscopes FritzscheJNCS70 ; FriedrichJAP00 ; MoralesJAP02 ;
LazarenkoAIP16 ; BurrJVST10 . We map the temperature contours of a nanowire
and an embedded micro-heater where the contour areas differ by three orders of
magnitude.

To demonstrate the versatility and practicality of this technique, a recording
head from a commercial hard disk drive is used. The head of the hard disk
drive provides a unique platform for such studies as it has several embedded
heat sources, which differ in heated area by three orders of magnitude WuAPL16
; ShimizuIEEE11 . At  the microscale, it has a micro-heater, which is used to
adjust the clearance between the head and the rotating disk RajauriaPRL18 .
The micro-heater is embedded a few micrometers from the surface, and  it
produces a microscale temperature contour. At  the nanoscale, it has a nano-
heater, which is used both as a heater and a thermometer. The nano-heater
consists of a 200 nm  wide, 1 $\mu$m long, and 20 nm thick metal wire that is
embedded 2 nm from the surface. Figure 1(c) shows the simulated surface
temperature contours for the micro-heater and  the nano-heater. In comparison
to the micro-heater ’s temperature field, the nano-heater’s areal temperature
map is three orders of magnitude smaller (see Figure 1(d)).  Note that the
nano-heater has a temperature coefficient of resistance (TCR) of 0.003/K,
which is used to measure the average surface temperature. In this paper, it is
quantitatively demonstrated that the temperature measured from the phase
change temperature contour (PCTC) technique agrees well with the measured
average surface temperature and the thermal simulation for both the micro-
heater and the nano-heater.

## Results and Discussions

### Self-heating of the nano-heater

Figure 2: Self-heating of the nano-heater. (a) Schematic diagram of the nano-
heater. Nanowire with dimension 1 $\mu$m $\times$ 20 nm $\times$ 200 nm is
electrically connected to two pads. (b)-(f)  Show the AFM  images of the
device at different micro-heater bias condition s. (c)-(f)  Show the
depression in the topography from the phase transition around the nano-heater.
(g) Shows the constructed temperature contour from the PCTC technique and the
simulation for the nano-heater power of 1.37 mW. (h) Shows the measured phase
change transition area as a function of dissipation power in the nano-heater
power. The red dash line corresponds to the simulation of an isotherm contour
for the glass transition temperature Tg. In the inset, the black dash line
shows the estimated average surface temperature along the nano-heater from the
resistance change in the nano-heater and the measured isotherm from the PCTC
technique (Red dot). Estimated error bar in average surface measurement is
0.04 oC.  Figure 3: Time response of $Ge_{2}Sb_{2}Te_{5}$ at a constant nano-
heater power of 0.68 mW. (a) The transition area of $Ge_{2}Sb_{2}Te_{5}$ phase
change with the the accumulated heating time. (b) Shows the temperature map of
the nano-heater. (c) Shows the temperature as a function of the hotspot width
across the nano-heater.

To demonstrate the technique, we first characterize the Joule heating of the
nano-heater inside the head of the hard disk drives. The head surface is
coated with a  22 nm thick layer of $Ge_{2}Sb_{2}Te_{5}$ thin film. The nano-
heater is biased using a current source,  across which the measured voltage
drop is used to estimate  the resistance increase of the nano-heater due to
the dissipated Joule heat (Figure 2(a)).  The resistance change of the nano-
heater is used to estimate the average temperature increase  of the device
temperature by $R_{T}=R_{0}(1+\alpha\Delta T)$, where  $R_{0}$ is the room-
temperature resistance at low current bias where no significant self-heating
occurs, $R_{T}$ is the resistance at the bias corresponding to the temperature
T, $\alpha$ is the temperature coefficient of resistance (TCR), and $\Delta$T
is the average temperature rise due to the Joule  heating. The temperature
coefficient of resistance ($\alpha$) 0.003/K is determined separately in an
oven using a 4-probe measurement scheme (see Supplementary section 1).  Note
that the effect of the thin layer on the heat transport of the system is
negligible (see Supplementary section 2).

The amorphous $Ge_{2}Sb_{2}Te_{5}$ is a chalcogenide phase change material
that crystallizes at T${}_{g}\sim 149^{\circ}$C  for a dwell time of 5
minutes. This crystallization is accompanied by an increase in density and
volume reduction,  where AFM topography measurement shows as a reduction in
the film height. Figures 2(b-f) show the AFM topography micrograph s
corresponding to different  powers in the nano-heater. For nano-heater power
smaller than 0.60 mW, the AFM shows no change in the topography of the
$Ge_{2}Sb_{2}Te_{5}$ film over the nano-heater, indicating that the surface
temperature is lower than the crystallization temperature.  When the nano-
heater power is 0.75 mW, a small depression in the topography is observed
centered at the hot spot of the nano-heater as shown in Figure 2(c). A further
increase in the nano-heater power leads to a gradual increase in the area
undergoing the crystallization, which is indicated by the lateral growth of
the depressed area in the AFM images.  Note that the boundary of the
topography depression corresponds to the isotherm of the crystallization
temperature. The evolution of the temperature contour area agrees reasonably
well with the simulation  as shown in Figure 2(h). Furthermore, we use the
transition boundary measured at different nano-heater powers to map the
temperature of the device. The rate of phase transition in
$Ge_{2}Sb_{2}Te_{5}$ is a function of both the temperature and the time. Here,
the power in the nano-heater is increased incrementally with a fixed dwell
time until the initial transition boundary is observed. The last transition
boundary corresponds to the calibration temperature Tg at the largest heater
power Po=1.37 mW (Figure 2(f)). Assuming that the temperature is linear with
the applied power, the temperature isotherm Ti at each previous transition
boundary (Figures 2(c-e)) is given by:

$T_{i}=T_{g}\frac{P_{o}}{P_{i}}$ (1)

where, Tg is the calibrated $Ge_{2}Sb_{2}Te_{5}$ crystalline transition
temperature for the dwell time of 300 seconds during which the nano-heater is
powered on, Po is the nano-heater power at which the final transition boundary
is measured, and Pi is the previous power with Pi $<$ Po in the nano-heater.
Figure 2(g) shows the constructed temperature map of the device along with the
simulation for Po=1.37 mW.

Figure 4: Self-heating of the micro-heater. (a)-(d) Show optical micrograph s
of the device at different micro-heater powers. (b)-(d) Show the reflectivity
increase at the center of the micrographs corresponding to the phase change
due to the temperature rise of the micro-heater. (e) Shows the average surface
temperature measured by the nano-heater, acting as a thermometer, along with
the  critical point for which the phase transition is measured from the
optical micrograph. Estimated error bar in average surface measurement is 0.04
oC. (f) Shows the constructed temperature map of the device along with the
simulation at the micro-heater power of 50 mW.

Figure 2(h) shows the phase change transition area calculated from the
topography depression in the AFM image as a function of the power dissipated
in the nano-heater. When the nano-heater power is lower than 0.60 mW, the
transition area is zero signifying that the surface temperature is lower than
the glass transition temperature Tg everywhere. At higher powers the
transition area grows linearly with the dissipated power in the nano-heater.
The simulation results are shown as the red dash line. Both the experiment and
the simulation show a sharp increase in the phase change transition area
beyond 0.60 mW. At much larger bias currents the poor match is because of our
inability to capture the exact structural details in the simulation  such as
the actual thermal boundary conditions and the various material parameters.
To further confirm the surface temperature, we simultaneously measure the
resistance of the nano-heater and  use a  TCR of 0.003/K to estimate the
average temperature of the heater. Figure 2(h) inset shows the measured
surface temperature of the nano-heater as a function of the power  dissipated.
The red dot is the temperature from the PCTC technique at 0.60 mW  derived
from the x-axis intercept of the Figure 2(h). As expected, the temperature
from the PCTC technique matches the temperature reading given by the
resistance change. The excellent agreement between the temperature measured
using the PCTC and TCR technique along with the simulation demonstrates that
the PCTC technique can precisely map the high operational temperature of the
nanoscale heater embedded in the chip.

Next we study the time response of the $Ge_{2}Sb_{2}Te_{5}$ film to construct
the temperature map of the device at a constant nano-heater power. Here, the
growth of the transition boundary is tracked over time. It should be noted
that the transient response of the nano-heater is six orders of magnitude
faster than the $Ge_{2}Sb_{2}Te_{5}$ phase change. Assuming that the phase
change conversion follows an Arrhenius model and the conversion is linear with
time, the temperature for each transition boundary is derived by monitoring
the time needed for each transition boundary to develop. The temperature
($T_{i}$) at time $t_{i}$ is given by:

$T_{i}(t_{i})=\left[-\frac{k_{B}}{E_{A}}\left(ln\frac{1}{t_{i}}-ln\frac{1}{t_{cal}}\right)+\frac{1}{T_{g}(t_{cal})}\right]^{-1}$
(2)

where, $k_{B}$ is the Boltzmann constant, $E_{A}$ $\sim$ 2.6 eV is the
activation energy $Ge_{2}Sb_{2}Te_{5}$ transition, T${}_{g}\sim$ 149∘C is the
calibrated crystallization temperature at dwell time tcal =300 sec. It is
worth noting that the temperatures of isotherms corresponding to shorter dwell
times ($t_{i}$ $<$ $t_{cal}$) are higher than Tg.

Figure 3(a) dots show the phase change transition area around the nano-heater
as a function of the accumulated time for a constant power of 0.68 mW (see
Supplementary section 3 and video TimeDependence.mp4). The red line is the
exponential fit with a time constant of 37.6 seconds. The transition
temperature at different accumulated heating time is determined by using Eq.
2. Figure 3(b) shows the constructed temperature map of the device at the
constant nano-heater power of 0.68 mW. In comparison to Figure 2(g), the
temperature map is smaller and more elliptical since the nano-heater power is
almost 50$\%$ smaller. Figure 3(c) shows the temperature across nano-heater as
a function of distance demonstrating the high resolution of the PCTC scheme.
The continuous nature of our thin film allows for a higher spatial resolution,
which is limited only by the grain size of $Ge_{2}Sb_{2}Te_{5}$ (sub 20 nm)
and the resolution of the imaging technique ($<$ 10 nm) BrintlingerNL08 .

### Self-heating of the micro-heater

To demonstrate the versatility of the PCTC technique, we now characterize the
Joule heating of a much larger micro-heater embedded in the head of the hard
disk drive.  The temperature contour of the micro-heater is three orders of
magnitude larger than  that of the nano-heater embedded in the same chip (as
shown in Figure 1)  in terms of the contour area.  The micro-heater is biased
using a current source and the measured voltage drop across the nano-heater
(thermometer) is used to estimate the dissipated Joule heat. The dwell time of
300 seconds at a constant micro-heater power is much longer than the thermal
response time of the heater and the phase transition time beyond which the
physical, optical and electrical properties change.  The temperature contour
of the micro-heater is mapped using an optical microscope by simply imaging
the reflectivity change in the transition area. Figures 4(a-d) show the
optical micrographs at different micro-heater power s. No change in
reflectivity is observed  below the dissipated power of 42 mW in the micro-
heater. At the power of 43 mW,  an increase in the reflectivity is observed at
the center of  the thermal hotspot due to the micro-heater.  Note that the
boundary of the transition area corresponds to the crystallization
temperature. In comparison to the nano-heater, here the micro-heater requires
60 times more power to achieve the same surface temperature since the micro-
heater is more deeply embedded and heats up a much larger volume. At higher
micro-heater power s, the growth in the transition area indicate s an increase
in  the thermal spot size with the same crystallization temperature Tg
isotherm.

Figure 4(f) shows the constructed temperature map of the device along with the
simulation.  The dimensions and the overall shape of the transition contour
from the experiments match well with the simulated Tg isotherms.  Furthermore,
we utilized the nano-heater as a thermometer by monitoring its resistance
change at a very low current of 0.1 mA, in order to avoid self-heating, to
measure the temperature rise due to the micro-heater. Figure 4(e) shows the
surface temperature measured by the nano-heater as a function of the power
dissipated in the micro-heater. The red line shows the micro-heater power
beyond which the phase transition is observed in the optical micrograph. The
expected rise of the surface temperature as derived from both the PCTC
technique and from the nano-heater (acting as ’thermometer’) is 2.9 K/mW. This
shows an excellent agreement between the PCTC technique, the measured surface
temperature and the simulation for the temperature map of the micro-heater.

Figure 5: Calibration. (a) Cartoon to show the phase transition in
$Ge_{2}Sb_{2}Te_{5}$. (b) Inset: Optical micrograph of the calibration sample.
Main: Thickness (measured using AFM, with vertical resolution of 0.05 nm) of
the calibration sample as a function of oven temperature. The glass transition
temperature Tg is 149∘C.

### Temperature calibration

Precise and relatively simple temperature calibration is a key advantage of
the PCTC method compared to other techniques that require extensive
temperature calibration. To calibrate the crystallization temperature Tg of
the phase change material $Ge_{2}Sb_{2}Te_{5}$, we rely on  a pronounced
structural property change at the phase change condition (Figure 5(a)). The
rate of this transformation from an amorphous state to the crystalline rock
salt structure is well characterized by an activation energy of about 2.6 eV
FriedrichJAP00 ; MoralesJAP02 . Due to this activation energy driven process
the crystallization temperature is dependent on the dwell time of the sample.
For this reason, the dwell time at a constant power condition is fixed at 5
minutes for both the nano-heater and the micro-heater. This dwell time is much
longer than the response time of the two heaters SchreckIEEE92 . Figure 5(b)
inset shows the calibration sample, which is  a photo-lithography defined 22
nm thick $Ge_{2}Sb_{2}Te_{5}$ layer on top of a silicon wafer.  The main
figure shows the film thickness as a function of oven temperature. The dwell
time at a constant oven temperature is 5 minutes after which the sample is
allowed to cool down to room temperature, and the film thickness is measured
using the AFM. At T = 149 ∘C, the film thickness reduces, indicating a phase
transition from the amorphous state to the crystalline state, as confirmed by
the X-ray diffraction pattern (see Supplementary section 4).

The uncertainty in the temperature  derived from this PCTC technique is
primarily due to the fact that the crystallization rate of the phase change
material does not have a large abrupt jump at a single temperature (as in a
first order phase transition). As a result, the full temperature history of
the sample, not just the last power used, can influence the size of the
observed contour (see Supplementary section 5). For both the nano-heater and
the micro-heater, the estimated temperature step is kept at 10 K, which leads
to a slight under-prediction of  the temperature by around 2 K at the 149 ∘C
crystallization condition.

It is noteworthy that in most other  high resolution temperature mapping
techniques,  the application is limited to the large scale devices. Here we
demonstrated the thermal measurements of the two heaters that produced thermal
contours with dimensions that differ by three orders of magnitude, while
maintaining sub 20 nm resolution in both cases. This technique is extremely
versatile and  does not require the use of expensive microscopes like STEM. At
the microscale, even an inexpensive optical microscope can be used to map the
temperature of the hotspots or heat sources in operating microelectronics
devices.

Finally, the limitations of the presented PCTC technique are discussed.
Although the technique is versatile and  can be used for nanoscale to
microscale spatial heat sources with minimal calibration challenges, it has
two main limitations. First, the  areal contour represents the isotherm at Tg,
and to extract the temperature gradient one needs to perform mapping for at
least two different power levels  or track the phase change with time.
Secondly, the technique is limited to  a device temperature higher than the
crystallization temperature of the deposited phase change material. In our
case, $Ge_{2}Sb_{2}Te_{5}$ crystallizes at temperature of 149 ∘C for the dwell
time of 5 minutes is used. Both embedded heaters are able to reach
temperatures higher than this Tg. For other systems where reaching a similar
temperature would be difficult or impossible, the issue could be overcome by
choosing a different composition or phase change materials with lower Tg
XiongIEEE16 ; ChenIEEE98 ; FuRSC15 ; AlkanSE12 .

## Conclusion

To summarize, we introduce a versatile phase-change-material-based temperature
mapping technique for operational microelectronic devices  that can spatially
resolve temperature from nanoscale to microscale dimensions. It can be used to
characterize surface temperatures with neglectable temperature interference
due to the deposited measurement film and with minimal calibration. A thorough
understanding of the heat dissipation in various nanoscale devices, such as
the aforementioned nano-heater, may lead to more efficient and powerful
integrated chips, and hence holds great economic value to the industry.

## Methods

### Experimental set-up

The microelectronic device is held on a metal fixture with electrical pins.
The components inside the device such as the heaters are powered by Keithley
2602 SYSTEM SourceMeter, which is controlled by a Python script. A 22 nm thin
film of Ge2Sb2Te5 is sputtered on the surface of the device. The topography
change or the reflectivity change of the thin film due to the heaters are
measured by Digital Instruments Dimension 3100 AFM or an optical microscope
respectively. In the temperature calibration of the Ge2Sb2Te5, a silicon wafer
with a photo-lithography defined 22 nm thick Ge2Sb2Te5 layer is heated in a
customized copper chamber, where the temperature is measured by a type-K
thermocouple. The thickness of the layer is also measured by the AFM.

### Simulation

Thermal simulations of the micro-heater and nano-heater devices in Figure 1
(c) and 1(d) were performed using finite element models in ANSYS Mechanical
APDL version 17.2. The nano-heater thermal simulation in Figure 2 (g) was
modeled using a finite element model in the ANSYS Workbench Thermo-electric
module version 17.2. The device surfaces were modeled using a convection
cooling boundary condition with a coefficient of 50 W/(m2·K).

##  I Acknowledgements

The authors would like to thank I. McFadyen for careful reading of the
manuscript. Q.C would like to thank Western Digital Corporation for financial
support for the internship.

##  II Author Contributions

The work was conceived by S.R. and E.S. The experimental setup were designed
and implemented by Q.C. and S.R. The experiments were performed and analyzed
by Q.C. under the supervision of S.R. and E.S. Thermal simulation was done by
R.S. Phase change material was deposited by J.R. The manuscript was written by
S.R. and Q.C with comments and input from all authors.

##  III Author information

Corresponding author<EMAIL_ADDRESS>
The authors declare they have no competing financial interests.

## References

  * (1) Haensch, W., Nowak, E. J., Dennard, R. H., Solomon, P. M., Bryant, A., Dokumaci, O. H., Kumar, A., Wang, X., Johnson, J. B. and Fischetti, M. V. Silicon CMOS devices beyond scaling.  IBM J. Res. Dev. 50, 339-361 (2006).
  * (2) Cavin, R., Zhirnov, V., Herr, D., Avila, A. and Hutchby, J. Research directions and challenges in nanoelectronics. Journal of Nanoparticle Research 8, 841-858 (2006).
  * (3) Pop, E. Energy dissipation and transport in nanoscale devices. Nano Research 3, 147-169 (2010).
  * (4) Dejene, F. K., Flipse, J., Bauer, G. E. W. and van Wees, B. J. Spin heat accumulation and spin-dependent temperatures in nanopillar spin valves. Nature Physics 9, 636-639 (2013).
  * (5) Jin, C. Y., Li, Z., Williams, R. S., Lee, K.-C. and Park, I. Localized temperature and chemical reaction control in nanoscale space by nanowire array. Nano Letters 11, 4818-4825 (2011).
  * (6) Brites, C.D.S. , Lima, P.P. , Silva, N.J., Millan, A., Amaral, V.S. , Palacio, F. and Carlos L.D. Thermometry at the nanoscale. Nanoscale 4, 4799-4829 (2012).
  * (7) Lee, J. and Kotov, N. A. Thermometer design at the nanoscale. Nanotoday 2, 48-51 (2007).
  * (8) Wilde, Y.D., Formanek, F., Carminati, R., Gralak, B., Lemoine, P.-A., Joulain, L., Mulet, J.-P., Chen, Y. and Greffet J.-J. Thermal radiation scanning tunnelling microscopy. Nature 444, 740-743 (2006).
  * (9) Deshpande, V. V., Hsieh, S., Bushmaker, A.W., Bockrath, M. and Cronin S.B. Spatially resolved temperature measurements of electrically heated carbon nanotubes. Physical Review Letters 102, 105501 (2009).
  * (10) Tessier, G., Bardoux, M., Boue, C., Filloy, C. and Fournier D. Back side thermal imaging of integrated circuits at high spatial resolution. Applied Physics Letters 90, 171112 (2007).
  * (11) Nonnenmacher, M. and Wickramasinghe H.K. Scanning probe microscopy of thermal conductivity and subsurface properties. Applied Physics Letters 61, 168 (1992).
  * (12) Shi, L., Plyasunov, S., Bachtold, A., McEuen, P.L. and Majumdar A. Scanning thermal microscopy of carbon nanotubes using batch-fabricated probes. Applied Physics Letters 77, 4295 (2000).
  * (13) Mecklenburg, M., Hubbard, W.A., White, E.R., Dhall, R., Cronin, S.B., Aloni, S. and Regan, B.C. Nanoscale temperature mapping in operating microelectronic devices. Science 347, 629-632 (2015).
  * (14) Brintlinger, T., Qi, Y., Baloch, K. H., Goldhaber-Gordon, D. and Cumings J. Electron Thermal Microscopy. Nano Letters 8, 582-585 (2008).
  * (15) Guo, H., Khan, M. I., Cheng, C., Fan, W., Dames, C., Wu J. and Minor A. M. Vanadium dioxide nanowire-based microthermometer for quantitative evaluation of electron beam heating. Nature Communications 5, 4986 (2014).
  * (16) Gao Y. and Bando Y. Carbon nanothermometer containing gallium. Nature 415, 599 (2002).
  * (17) Martinek, J., Klapetek, P. and Campbell, A. C. Methods for topography artifacts compensation in scanning thermal microscopy. Ultramicroscopy 155, 55-61 (2015).
  * (18) Menges, F., Mensch, P., Schmid, H., Riel, H., Stemmer A. amd Gotsmann B. Temperature mapping of operating nanoscale devices by scanning probe thermometry. Nature Communications 7, 10874 (2016).
  * (19) Kinkhabwala, A., Staffaroni, M., Suzer, O., Burgos, S. and Stipe B. Nanoscale thermal mapping of HAMR heads using polymer imprint thermal mapping. IEEE Transactions on Magnetics 52, 1-4 (2015).
  * (20) Kim, K., Jeong, W., Lee, W. and Reddy, P. Ultra-high vacuum scanning thermal microscopy for nanometer resolution quantitative thermometry. ACS Nano 81, 114901 (2012).
  * (21) Hwang, G., Chung, J. and Kwon, O. Enabling low-noise null-point scanning thermal microscopy by the optimization of scanning thermal microscope probe through a rigorous theory of quantitative measurement. Rev. Sci. Instrum. 85, 114901 (2014).
  * (22) Okabe, K., Inada, N., Gota, C., Harada, Y., Funatsu T. and Uchiyama S. Intracellular temperature mapping with a fluorescent polymeric thermometer and fluorescence lifetime imaging microscopy. Nature Communications 3, 705 (2012).
  * (23) Kucsko, G., Maurer, P. C., Yao, N. Y., Kubo, M., Noh, H. J., Lo, P. K., Park, H. and Lukin M. D. Nanometre-scale thermometry in a living cell. Nature 500, 54-58 (2013).
  * (24) Idrobo, J. C., Lupini, A. R., Feng, T., Unocic, R. R., Walden, F. S., Gardiner, D. S., Lovejoy, T. C., Dellby, N., Pantelides, S. T. and Krivanek, O. L. Temperature Measurement by a Nanoscale Electron Probe Using Energy Gain and Loss Spectroscopy. Physical Review Letters 120, 095901 (2018).
  * (25) Feng, et al. Electronic Structure and Enhanced Charge-Density Wave Order of Monolayer VSe2. Nano Letters 18, 4493–4499 (2018).
  * (26) Christofferson, J., Maize, K., Ezzahri, Y., Shabani, J., Wang, X. and Shakouri A. Microscale and nanoscale thermal characterization techniques. J. Electron. Packag. 130, 041101 (2008).
  * (27) Menges, F., Riel, H., Stemmer, A. and Gotsmann, B. Quantitative thermometry of nanoscale hot spots. Nano Letters 12, 596-601 (2012).
  * (28) Fritzsche, H. and Ovshinsky S.R. Calorimetric and dilatometric studies on chalcogenide alloy glasses. Journal of Non-Crystalline Solids 2, 148-154 (1970).
  * (29) Friedrich, I., Weidenhof, V., Njoroge, W., Franz, P. and Wuttig M. Structural transformations of $Ge_{2}Sb_{2}Te_{5}$ films studied by electrical resistance measurements. Journal of Applied Physics 87, 4130 (2000).
  * (30) Morales-Sanchez, E., Prokhorov, E. F., Mendoza-Galvan, A. and Gonzalez-Hernandez J. Determination of the glass transition and nucleation temperatures in $Ge_{2}Sb_{2}Te_{5}$ sputtered films Journal of Applied Physics 91, 697 (2002).
  * (31) Lazarenko, P. I., Sherchenkov, A. A., Kozyukhin, S. A., Babich, A. V., Timoshenkov, S. P., Gromov, D. G., Shuliatyev, A. S. and Redichev E. N. Electrical properties of the $Ge_{2}Sb_{2}Te_{5}$ thin films for phase change memory application. AIP Conference Proceedings 1727, 020013 (2016).
  * (32) Burr G.W. Phase change memory technology. Journal of Vacuum Science and Technology B 28, 223 (2010).
  * (33) Wu, H., Xiong, S., Canchi, S., Schreck, E. and Bogy, D. Nanoscale heat transfer in the head-disk interface for heat assisted magnetic recording. Applied Physics Letters 108, 093106 (2016).
  * (34) Shimizu, Y., Xu, J., Kohira, H., Kurita, M., Shiramatsu, T.and M. Furukawa M. Nano-scale defect mapping on a magnetic disk surface using a contact Sensor. IEEE Transactions on Magnetics 35, 3426-3432 (2011).
  * (35) Rajauria, S., Ruiz, O., Canchi, S.V., Schreck, E. and Dai Q. Electrostatically tunable adhesion in a high speed sliding interface. Physical Review Letters 120, 026101 (2018).
  * (36) Schreck, E., Fontana, R. E., and Singh G.P. Thin film thermocouple sensors for measurement of contact temperatures during slider asperity interaction on magnetic recording disks. IEEE Transactions on Magnetics 28, 2548-2550 (1992).
  * (37) Xiong, F., Yalon, E., Behnam, A., Neumann, C.M., Grosse, K.L., Deshmukh, S. and Pop E. Towards ultimate scaling limits of phase-change memory. In Proc. IEEE International Electron Devices Meeting 4.1.14.1.4 (IEEE, San Francisco, CA, USA, 2016).
  * (38) Chen, Y. M. and Kuo P. C. Effect of Ag or Cu doping on erasable phase-change Sb-Te thin films. IEEE Transactions on Magnetics 34, 432-434 (1998).
  * (39) Fu, X., Kong, W., Zhang, Y., Jiang, L., Wangb J. and Lei, J. Novel solid–solid phase change materials with biodegradable trihydroxy surfactants for thermal energy storage. RSC Adv. 5, 68881-68889 (2015 ).
  * (40) Alkan, C., Gunther, E., Hiebler, S., Ensari O. F. and Kahraman D. Polyurethanes as solid–solid phase change materials for thermal energy storage. Solar Energy 86, 1761-1769 (2012).