# The ALMA View of Positive Black Hole Feedback in the Dwarf Galaxy Henize
2–10

Hansung B. Gim Department of Physics, Montana State University, Bozeman, MT
59717, USA  Amy E. Reines Department of Physics, Montana State University,
Bozeman, MT 59717, USA  eXtreme Gravity Institute, Department of Physics,
Montana State University, Bozeman, MT 59717, USA Hansung B. Gim
<EMAIL_ADDRESS>
###### Abstract

Henize 2–10 is a dwarf starburst galaxy hosting a $\sim 10^{6}~{}M_{\odot}$
black hole (BH) that is driving an ionized outflow and triggering star
formation within the central $\sim 100$ pc of the galaxy. Here we present ALMA
continuum observations from 99 to 340 GHz, as well as spectral line
observations of the molecules CO (1–0, 3–2), HCN (1–0, 3–2), and HCO$+$ (1–0,
3–2), with a focus on the BH and its vicinity. Incorporating cm-wave radio
measurements from the literature, we show that the spectral energy
distribution of the BH is dominated by synchrotron emission from 1.4 to 340
GHz with a spectral index of $\alpha\approx-0.5$. We analyze the spectral line
data and identify an elongated molecular gas structure around the BH with a
velocity distinct from the surrounding regions. The physical extent of this
molecular gas structure is $\approx 130~{}{\rm pc}\times 30$ pc and the
molecular gas mass is $\sim 10^{6}~{}M_{\odot}$. Despite an abundance of
molecular gas in this general region, the position of the BH is significantly
offset from the peak intensity, which may explain why the BH is radiating at a
very low Eddington ratio. Our analysis of the spatially-resolved line ratio
between CO J=3–2 and J=1–0 implies that the CO gas in the vicinity of the BH
is highly excited, particularly at the interface between the BH outflow and
the regions of triggered star formation. This suggests that the cold molecular
gas is being shocked by the bipolar outflow from the BH, supporting the case
for positive BH feedback.

galaxies:active – galaxies: dwarf –galaxies: individual (He 2–10) - ISM:
molecules

††facilities: ALMA, VLA, LBA, HST††software: CASA (Bean et al., 2022), AIPS
(Greisen, 2003), CARTA (Comrie et al., 2021), R (R Core Team, 2022), 3D-Barolo
(Di Teodoro & Fraternali, 2015)

OIII

##  1 Introduction

Galaxy evolution is governed by the balance between star formation and
quenching (Somerville & Davé, 2015), and Active Galactic Nuclei (AGN) are
thought to play a crucial role in the extinguishing of star formation in
massive galaxies. Prior to the last decade, massive black holes (BHs) that
sometimes accrete and shine as AGNs were almost exclusively found in giant
galaxies (e.g., Kormendy & Ho, 2013). However, we now know that AGNs in dwarf
galaxies are much more common than previously thought and the BHs have typical
masses of $M_{\rm BH}<10^{6}~{}M_{\odot}$ (for reviews, see Greene et al.,
2020; Reines, 2022). Searching for and studying these systems is important for
our understanding of the role of AGN fueling and feedback in the low-mass
regime and may also provide clues about the origin of the first BH “seeds” in
the early Universe (e.g., Volonteri et al., 2021).

Henize 2–10 (He 2–10) is a nearby (d$=$9 Mpc) dwarf starburst galaxy (e.g.,
Nguyen et al., 2014) that is particularly remarkable for displaying evidence
of positive AGN feedback (i.e., enhancing star formation; Schutte & Reines,
2022, see Figure 1 here). Multiwavelength studies have shown that this galaxy
hosts a massive BH in its center (Reines et al., 2011) through coincident
peaks at radio and X-ray wavelengths (Reines et al., 2016) detected with the
Very Large Array (VLA) and the Chandra X-ray Observatory (Chandra),
respectively. The BH is also detected as a strong non-thermal source by high-
angular resolution observations using the Long Baseline Array (LBA, Reines &
Deller, 2012). The mass of the central BH is estimated to be $\rm M_{BH}\sim
10^{6}\;M_{\odot}$ based on stellar velocity measurements (Riffel, 2020) and
the scaling relation between total stellar mass and BH mass (Reines &
Volonteri, 2015). Combined with the X-ray luminosity, this indicates the BH is
accreting at a very low Eddington ratio (Reines et al., 2016). A recent study
using Hubble Space Telescope (HST) spectroscopy revealed that the central BH
is driving an outflow of ionized gas, which is triggering the formation of
star clusters located $\sim$70 pc ($\sim$ 1$\farcs$5) from the BH (Schutte &
Reines, 2022). This is the first case of positive BH feedback found in a dwarf
galaxy and may represent a low-power analog of radio-loud AGNs experiencing
“jet-mode” feedback (Schutte & Reines, 2022).

Figure 1: The four color composite image of He 2–10 is generated from F330W
(B), F606W (V), F658N (Hα), and F814W (I) from the HST. The white box shows
the central 6″$\times$4″ region, and the inset shows the H$\alpha$$+$continuum
image presenting the outflow from the BH. (Credit to NASA, ESA, Zachary
Schutte, & Amy Reines)

Molecular gas and dust can both feed and obscure AGNs. While the cold gas
distributions of He 2–10 have been extensively studied for the galaxy as a
whole and with a focus on the star-forming regions, very little is known about
the molecular gas in the immediate vicinity of the BH. Here we present and
analyze high-resolution Atacama Large Millimeter/submillimeter Array
(ALMA)111The National Radio Astronomy Observatory is a facility of the
National Science Foundation operated under cooperative agreement by Associated
Universities, Inc. observations of He 2–10 with a specific focus on the
central massive BH and its immediate surroundings ($\sim$ 130 pc). A summary
of previous studies of the molecular gas in He 2–10 is given below.

Baas et al. (1994) first observed the three CO transitions (J=1–0, J=2–1, and
J=3–2) of He 2–10, showing that hot and cold molecular gas is distributed
throughout the galaxy. The first interferometric observations of CO J=1–0 with
a synthesized beam size of 6$\farcs$5$\times$5$\farcs$5 (284 pc$\times$240 pc)
revealed an elongated morphology, indicating that He 2–10 could be a post-
merger system and/or experiencing tidal interactions (Kobulnicky et al.,
1995). Meier et al. (2001) observed the CO J=3–2 line to investigate the
interstellar medium (ISM) conditions in the central 13″ (567 pc) and found a
gas temperature of $\rm T_{K}\sim 5-10$ K and a number density of $\rm
n(H_{2})>10^{3.5}$ cm-3. Higher transitions of CO lines, J=3–2, 4–3, 6–5, and
7–6, were observed by Bayet et al. (2004). They found that the central region
hosts warm and dense molecular gas with T${}_{K}\approx 50-100$ K and $\rm
n(H_{2})>10^{4}$ cm-3, which is an order of magnitude higher than the
temperature calculated with lower transitions of CO lines (J = 1–0, 2–1 and
3–2) by Meier et al. (2001). Beck et al. (2018) observed the CO J=3–2 line
with a synthesized beam size of 0$\farcs$3 (13 pc) using ALMA and found
several clumpy filaments separated spatially and kinematically. Their
kinematic analysis led to the conclusion that the molecular filaments feed the
star formation in the star clusters. Imara & Faesi (2019) traced the CO J=1–0
line with a synthesized beam size of 0$\farcs$67$\times$0$\farcs$58 (29
pc$\times$25 pc) using ALMA and showed that He 2–10 has a more complex
morphology with clumps, ring-shaped structures, and high-density peak
emissions.

Higher critical density molecular gas tracers have also been observed in He
2–10. Imanishi et al. (2007) observed the HCN J=1–0 and HCO+ J = 1–0 line with
a beam size of 10$\farcs$8$\times$5$\farcs$5 (471 pc$\times$240 pc) and
reported the ratio of HCN (1–0)/HCO+ (1–0) $<0.6$, which is low compared to
Luminous InfraRed Galaxies (LIRGs). The strong HCN J=1–0 line was also
recovered by the single dish observations of Santangelo et al. (2009). Johnson
et al. (2018) used ALMA observations tracing HCN J=1–0, HCO+ J=1–0, HNC J=1–0,
and CCH J=1–0 with a synthesized beam size of 1$\farcs$7$\times$1$\farcs$6 (74
pc$\times$70 pc) and showed that HCO+ J=1–0 correlates with the thermal radio
continuum in the intense central star forming regions.

Here we present the first study of the molecular gas in the immediate vicinity
of the central BH using high-resolution observations with ALMA. We describe
the observations and data reduction in § 2, the analysis and results in § 3,
and our conclusions in § 4.

##  2 ALMA Observations and Data Reduction

He 2–10 has been observed by ALMA in different bands for six different
projects. The specifics of the ALMA observations are summarized in Table 1,
encompassing essential information such as project codes, bands, frequency
ranges, spectral lines, total integration time, bandpass, phase, and amplitude
calibrators. The names of PIs are presented in the caption located below the
table. We downloaded all the data listed in Table 1 from the ALMA archive
website222The ALMA archival website is at https://almascience.nrao.edu/aq/.
The calibrated data were processed using the scriptForPI.py with the
appropriate version of the Common Astronomy Software Application (CASA, Bean
et al., 2022) as outlined in the README files. We also flagged some bad
antennas, which improved the sensitivity.

Table 1: ALMA observations of He 2–10

(1) | (2) | (3) | (4) | (5) | (6) | (7) | (8)
---|---|---|---|---|---|---|---
Project code | Band | Frequency | Lines | Integration | Bandpassa | Phaseb | Amplitudec
|  | (GHz) |  |  |  |  |
2011.0.00348.S (1) | 3 | 86.95 – 102.47 | HCN (1–0) | 204m 37s | J053851-440507 | J082601-223027 | Titan
HCO+ (1–0) | Ceres
2012.1.00413.S (2) | 6 | 248.79 – 267.73 | HCN (3–2) | 286m 47s | J0750+1231 | J0747-3310 | J1037-295
HCO+ (3–2) | J0538-4405 | J0826-2230 | J0538-4405
J1058+0133
2015.1.01569.S (3) | 3 | 111.99 – 115.18 | CO (1–0) | 187m 50s | J1037-2934 | J0826-2230 | J1107-4449
2016.1.00027.S (3) | 3 | 99.72 – 115.17 | J0538-4405 | J0846-2607
2016.1.00492.S (4) | 7 | 332.06 – 348.00 | CO (3–2) | 50m 29s | J0538-4405 | J0846-2607 | J0538-4405
2019.1.01641.S (5) | 3 | 96.91 – 112.51 | – | 103m 43s | J0725-0054 | J0846-2607 | J0725-0054
J1037-2934 | J1037-2934

Note. — $(a)$ Bandpass calibrators ; $(b)$ Phase calibrators ; $(c)$ Amplitude
calibrators

PI names: (1) Kelsey Johnson ; (2) Amy Reines ; (3) Nia Imara ; (4) Sara Beck
; (5) Zhiyu Zhang

###  2.1 Spectral Line Image Cubes

Spectral line image cubes were generated using CASA 6.5.1 with the continuum-
subtracted spectral line visibilities. The continuum subtraction was performed
using the CASA task uvcontsub, which uses line-free channels and a linear fit
of order 1 to estimate the continuum. The continuum-subtracted visibilities
were then imaged using the CASA task tclean with the hogbom algorithm for
deconvolution of the point spread function (PSF), briggsbwtaper, and
perchanweightdensity. These options provide a flat root mean square (RMS)
noise and beam size distribution over the channels, and the beam size
comparable to that of the continuum image. The velocity frame is set to BARY,
which normalizes the velocity to the barycentric system. The cleaning process
was repeated until the maximum of the residual map reached 2$\sigma$, where
$\sigma$ is the RMS noise of the image cube per channel. The line-free
channels for the continuum subtraction, image size, cell size, velocity
spacing, robust value for the briggsbwtaper, and rest frequency are summarized
in Table 2. The resulting image cubes were smoothed to the largest beam size
within the selected velocity range including lines with the CASA task
imsmooth. The resulting beam sizes and RMS noises of molecular line image
cubes are summarized in Table 2.

Table 2: Parameters for Imaging and Image Statistics

(1) | (2) | (3) | (4) | (5) | (6) | (7) | (8) | (9)
---|---|---|---|---|---|---|---|---
Line | Line-free | Image size | Cell size | Velocity spacing | robust | Rest frequency | Synthesized | RMS
Channels | beam size
|  | (pixel$\times$pixel) | (″ pixel-1) | (km s-1) |  | (GHz) |  | ($\mu$Jy beam-1)
HCN (1–0) | 301–750 | 512$\times$512 | 0.25 | 10.0 | 0.5 | 88.6316 | 1$\farcs$78$\times$1$\farcs$68 | 350a
1150–3539
HCO+ (1–0) | 300–3300 | 512$\times$512 | 0.25 | 10.0 | 0.5 | 89.1885 | 1$\farcs$84$\times$1$\farcs$68 | 350a
3500–3539
HCN (3–2) | 200–1600 | 1024$\times$1024 | 0.025 | 20.0 | 1.5 | 265.886 | 0$\farcs$19$\times$0$\farcs$16 | $<$ 470a
2000–3640
HCO+ (3–2) | 200–1600 | 1024$\times$1024 | 0.025 | 7.5 | 1.5 | 267.558 | 0$\farcs$18$\times$0$\farcs$16 | 325a
2500–3640
CO (1–0) | 100–1700 | 1024$\times$1024 | 0.05 | 1.693 | 0.0 | 115.208201 | 0$\farcs$45$\times$0$\farcs$28 | 1030a
2150–3740
CO (3–2) | 10–90 | 1024$\times$1024 | 0.05 | 1.693 | 0.5 | 345.79599 | 0$\farcs$31$\times$0$\farcs$28 | 1000a
161–230
99 GHz | – | 1024$\times$1024 | 0.05 | – | 0.0 | – | 0$\farcs$55$\times$0$\farcs$32 | 13.0
113 GHz | – | 1024$\times$1024 | 0.05 | – | 0.0 | – | 0$\farcs$45$\times$0$\farcs$27 | 18.0
251 GHz | – | 1024$\times$1024 | 0.03 | – | 1.5 | – | 0$\farcs$18$\times$0$\farcs$15 | 12.4
340 GHz | – | 1024$\times$1024 | 0.04 | – | 0.5 | – | 0$\farcs$32$\times$0$\farcs$27 | 27.7

  * •

a The unit is ($\mu$Jy beam-1 (channel width)-1) for spectral image cube.

###  2.2 Continuum Imaging

We created four different continuum images: low-frequency Band 3 (low-Band 3),
high-frequency Band 3 (high-Band 3), Band 6, and Band 7 continuum images. The
Band 3 continuum data are from ALMA project codes 2011.0.00348.S,
2015.1.01569.S, 2016.1.00027.S, and 2019.1.01641.S, which were split into low-
and high-frequency data at 110 GHz. This frequency was chosen to maximize the
angular resolution of the high-frequency Band 3 continuum image. The Band 6
and 7 continuum data are from ALMA project codes 2012.1.00413.S and
2016.1.00492.S, respectively.

The continuum images were created by combining the visibilities from the
continuum spectral windows (SPWs) with wide channel widths of 15.625 MHz and
the line-free channels from the SPWs containing spectral lines with narrow
channel widths (as explained in Section 2.1). The latter were smoothed to a
channel width of 15.625 MHz in SPWs for the continuum using the CASA task
mstransform. The continuum images were generated using the CASA 6.5.1 task
tclean with the briggs weighting function and the hogbom deconvolution
algorithm. The standard multi-frequency synthesis (MFS) was applied for
imaging the continuum data, except for the low-Band 3. For the low-Band 3
data, we applied the multi-term MFS (Rau & Cornwell, 2011), as its fractional
bandwidth is about 23.7% (with a bandwidth of 23.4 GHz, or 87.1–110.5 GHz,
relative to the effective frequency of 98.8 GHz). We used nterms=3, which
yields a more accurate intensity map by calculating the mean intensity,
spectral index, and curvature over the frequency. The cleaning was stopped
when it reached 2$\sigma$, where $\sigma$ is the RMS noise of the image. We
summarized the parameters used in the CASA task tclean and the resulting
synthesized beam sizes and RMS noise in Table 2.

The synthesized beam sizes and sensitivities of the continuum images in this
paper are slightly different from those obtained in Costa et al. (2021) for
113 and 340 GHz, while the characteristics of the 251 GHz image are the same.
Our 113 GHz continuum image was obtained by concatenating the data of
2015.1.01569.S, 2016.1.00027.S, and 2019.1.01641.S, which is different from
Costa et al. (2021) who used 2016.1.00027.S for the 113 GHz continuum image.
The synthesized beam size and RMS noise were 0$\farcs$38$\times$0$\farcs$21
and 17 $\mu$Jy beam-1, respectively (Costa et al., 2021), where the beam size
is 18% better than our 113 GHz continuum image. This difference is attributed
to adding two other observations to our data, 2015.1.01569.S and
2019.1.01641.S, since the robust value used in CASA task tclean was the same
as ours. The project codes 2015.1.01569.S and 2019.1.01641.S have maximum
baselines of 1.4 km and 181.9 m, respectively, while the project
2016.1.00027.S has a maximum baseline of 2.6 km. Adding these two data with
shorter baselines degraded the angular resolution when the data of
2016.1.00027.S is used solely. While the synthesized beam sizes are almost
same for the 340 GHz images, the RMS noises are 27.7 and 52.0 $\mu$Jy beam-1
for our image and that of Costa et al. (2021). The reason for this difference
is unclear.

###  2.3 Astrometry

It is important to facilitate the comparison of multi-wavelength observations
of Henize 2–10 taken with different telescopes. The location of the BH in He
2–10 was originally determined in the reference frame of 2MASS, due to its
absolute positional accuracy of $\lesssim$ 0$\farcs$1\. Images from HST, the
VLA, and Chandra were shifted to match the locations with the astrometric
frame of it 2MASS (Reines et al., 2011). Here, we identified the offsets in
our ALMA observations with respect to the 2MASS reference frame and applied
them to the ALMA observations.

We used the highest angular resolution ALMA continuum image at 251 GHz (with a
synthesized beam size of 0$\farcs$18$\times$0$\farcs$16) and the 1.4 GHz
continuum image observed at the Long Baseline Array (LBA) with a synthesized
beam size of 0$\farcs$11$\times$0$\farcs$03 to estimate the offset. The
astrometry of the LBA 1.4 GHz continuum image was corrected to that of the
2MASS (Reines & Deller, 2012). The shift was determined by matching the local
maximum in the 251 GHz continuum image to the center of the BH in the 1.4 GHz
continuum image. We estimated the shift of $\rm\Delta RA=$0$\farcs$0735 West
and $\rm\Delta DEC=$0$\farcs$095 North and applied this shift to all ALMA
spectral line image cubes and continuum images.

##  3 Analysis and Results

###  3.1 Continuum Emission and Radio SED of the Black Hole

Figure 2: Continuum images at 99 GHz (A), 113 GHz (B), 251 GHz (C), and 340
GHz (D). The central BH is marked with the black ellipse and red arrow in the
99 GHz continuum image (panel A). White contours represent the 99 GHz
continuum at levels of 50, 100, 150, and 200 $\mu$Jy beam-1. The beam size in
each image is shown as the gray ellipse at the bottom right.

Figure 2 presents the sub-mm continuum maps at 99 (first row), 113 (second
row), 251 (third row), and 340 GHz (fourth row), respectively. The white
contours show the 99 GHz continuum at 100, 150, 180, 200, and 250 $\mu$Jy
beam-1. The local maxima of the contours in the center indicate the position
of BH, which coincides with the positions determined by the VLA & LBA (Reines
et al., 2011; Reines & Deller, 2012). The continuum images show two prominent
blobs and one central source. The two strong emission blobs correspond to the
starburst regions, and the central source indicates the massive BH (Reines et
al., 2011; Reines & Deller, 2012). The position of the BH is defined in the
continuum images as the local maxima in the central region between the two
prominent starburst regions, as indicated in the top panel of Figure 2. The
two starburst regions are clearly visible in all the continuum maps; however,
the central BH is not seen in the 340 GHz continuum map. In the 340 GHz
continuum image, the BH is located between two small continuum sources located
0$\farcs$217 (10 pc) and 0$\farcs$526 (23 pc) away.

An earlier study by Allen et al. (1976) showed that He 2–10 as a whole is
dominated by nonthermal emission. However, high angular resolution VLA and
ALMA observations have shown that the central few hundred parsecs are
dominated by free-free emission from bursty star formation (Johnson &
Kobulnicky, 2003; Johnson et al., 2018; Costa et al., 2021). The central BH,
on the other hand, is governed by nonthermal synchrotron emission (Johnson &
Kobulnicky, 2003; Reines et al., 2011; Reines & Deller, 2012; Costa et al.,
2021). Although the VLA + ALMA SED analysis for the BH has been done by Costa
et al. (2021), we revisited this analysis with revised flux density
measurements of the submillimeter (sub-mm) continuum and inclusion of 1.4 GHz
flux density observed at the LBA (Reines & Deller, 2012).

The radio SED of the BH from 1.4 to 340 GHz is shown in Figure 3. In addition
to our new ALMA measurements in the sub-mm regime (see below), we include
measurements from the literature including 1.4 GHz data from the LBA (Reines &
Deller, 2012), 5 and 8 GHz data from the VLA (Reines et al., 2011), as well as
15, 22, and 33 GHz data from the VLA (Costa et al., 2021). Here, the flux
density at 33 GHz is not included in this analysis because the data suffered
from decorrelation due to weather conditions and the low declination of He
2–10 at this high frequency (Costa et al., 2021). Therefore, the flux density
at 33 GHz is presented with a lower limit (upper arrow) in Figure 3.

The flux densities of the BH in the sub-mm regime were estimated with our
reimaged ALMA continuum data. First, to help mitigate the impact of different
beam sizes, we convolved all the continuum images to have the same synthesized
beam size as the lowest one (0$\farcs$55$\times$0$\farcs$32 for the 99 GHz
image) using the CASA task imsmooth. While the inherent limitation of
interferometry may still result in different spatial sensitivities among the
observations, as well as the omission of large-scale emission, this would not
be a major concern if the emission associated with the BH is point-like.
However, if the emission is marginally resolved, we may be missing flux
density that could, for example, lead to the low measurement at 251 GHz seen
in Figure 3.

We performed aperture photometry on the sub-mm continuum using an elliptical
aperture equivalent to one beam size centered on the BH position (i.e.,
0$\farcs$55$\times$0$\farcs$32 with PA$=$-64.0∘ marked with a blue ellipse in
Figure 3). The arrow points to the BH which is shown as a black ellipse. The
measured flux densities of the BH are 140$\pm$14, 134$\pm$21, 29$\pm$46, and
57$\pm$48 $\mu$Jy at 99, 113, 251, and 340 GHz, respectively. Here, it is
worth noting that the uncertainty of the flux density may be underestimated
because we measured it as a mean of the RMS noise within an aperture. The RMS
noise maps were produced with the AIPS task RMSD (Greisen, 2003), which
calculates the RMS of signals after rejecting signals greater than 3$\sigma$
within 150 pixels of each pixel.

Figure 3: Continuum images at 99 GHz (A), 113 GHz (B), 251 GHz (C), 340 GHz
(D), and the SED analysis (E). Continuum images were convolved to have the
same synthesized beam size of 0$\farcs$55$\times$0$\farcs$32\. The central BH
is marked with the small black ellipse indicated by the arrow. The blue
ellipse is the aperture of our flux density measurements. The beam size in
each image is shown as the gray ellipse at the bottom right. The source A and
B in the 340 GHz continuum image (panel D) are adjacent sources of the BH. In
panel (E), the flux densities are depicted as black data points, with a black
arrow indicating the lower limit at 33 GHz. The best-fit model is represented
by a red line, and the 1$\sigma$ range of the best-fit model is delineated by
the red shaded region, providing a measure of the model’s uncertainty.

Our measurements are consistent with those of Costa et al. (2021) at 113 GHz,
but not at 251 and 340 GHz. Costa et al. (2021) reported flux densities of
70$\pm$30, 130$\pm$30 $\mu$Jy at 258 and 340 GHz, respectively, while our
measurements are 29$\pm$46 and 57$\pm$48 $\mu$Jy at these frequencies. At 251
GHz, our photometry shows 2.4 times less flux density than that of Costa et
al. (2021). Costa et al. (2021) likely measured the flux density of a point
source shown within the contour in the center of panel (C) of Figure 2.
However, this point source is smeared after convolution to the beamsize of
0$\farcs$55$\times$0$\farcs$32 in our measurement, as shown in panel (C) of
Figure 3. Since we applied aperture photometry on the convolved image, we
measured a smaller flux density. At 340 GHz, our measurement is 2.3 times
fainter than the flux density obtained by Costa et al. (2021). This difference
is due to the difference in measurements. There are two continuum sources near
the BH position with distances of 0$\farcs$22 (source A) and 0$\farcs$53
(source B) from the BH, as shown in panel (D) of Figure 3. Their flux
densities were measured as 125$\pm$80 and 342$\pm$130 $\mu$Jy for the sources
A and B, respectively, using PyBDSF (Mohan & Rafferty, 2015). Costa et al.
(2021) likely reported the flux density of source A as that of the BH. Our
aperture for the BH does not include source A and source B entirely, but
includes a part of source A, as seen in panel (D) of Figure 3. Therefore, our
reported values are more robust than the ones reported by Costa et al. (2021).

Our measurements have larger uncertainties than those of Costa et al. (2021)
at both 251 and 340 GHz. These larger uncertainties are due to the convolution
to the beamsize of 0$\farcs$55$\times$0$\farcs$32\. The differences in beam
area between the original and convolved images are 6.5 and 2.0 at 251 and 340
GHz, respectively. Convolution to the larger beam size introduces a larger
uncertainty. Quantitatively, the 251 GHz image with the original synthesized
beam size of 0$\farcs$18$\times$0$\farcs$14 has an RMS noise of 12.4 $\mu$Jy
beam-1, but its RMS noise increases to 24.7 $\mu$Jy beam-1 after convolution
to a synthesized beamsize of 0$\farcs$55$\times$0$\farcs$32\. Similarly, the
340 GHz image has an RMS noise of 27.7 $\mu$Jy beam-1 originally, but
increases to 38.7 $\mu$Jy beam-1 after convolution. Our measurements at 251
and 340 GHz had larger uncertainties, but it is important to match the
beamsizes across the frequencies in order to accurately measure the spectral
index. Additionally, it is important to measure the flux density of the BH at
the correct position, rather than the flux density of a nearby yet different
source. As a result, the flux densities we report are more accurately measured
than previous values.

In order to estimate the best-fit parameters and corresponding uncertainties
for an SED ranging from 1.4 to 340 GHz, we employed the R function nls for the
non-linear least square fitting and utilized confint2 to estimate the
uncertainties of the best-fit parameters. Four models were tested, including
i) a three-component model with non-thermal, free-free, and dust emission,
$S\sim S_{0,nth}\nu^{\alpha_{nth}}+S_{0,ff}\nu^{-0.1}+S_{0,dust}\nu^{4}$, ii)
a two-component model with non-thermal and free-free emission, $S\sim
S_{0,nth}\nu^{\alpha_{nth}}+S_{0,ff}\nu^{-0.1}$, iii) a two-component model
with non-thermal and dust emission, $S\sim
S_{0,nth}\nu^{\alpha_{nth}}+S_{0,dust}\nu^{4}$, and iv) a one-component model
with non-thermal emission, $S\sim S_{0,nth}\nu^{\alpha_{nth}}$. The models
that included dust emission did not yield physical results. The two-component
model with non-thermal and free-free emission had large uncertainties and
poorly constrained solutions. The one-component model with non-thermal
emission provided a meaningful solution, with the best-fit model of
$S=(1349\pm 217)\nu^{-0.52\pm 0.06}$ $\mu$Jy, represented by the red line with
a light red shaded region in panel E of Figure 3.

###  3.2 Molecular Line Emission

####  3.2.1 Moment Maps

Figure 4: Moment maps of CO (1–0), CO (3–2), HCN (1–0), HCO+ (1–0), and HCO+
(3–2) from top to bottom. Total intensity, velocity field, and velocity
dispersion maps are shown in the left, middle, and right panels, respectively.
The black contours indicate the flux density of the 99 GHz continuum,
corresponding to 50 (5$\sigma$), 100 (10$\sigma$), 150 (15$\sigma$), and 200
(20$\sigma$) $\mu$Jy per beam size of 0$\farcs$55$\times$0$\farcs$3\. The BH
is located in the center of the innermost black contour of the continuum. The
beam sizes are marked with solid ellipses in the bottom right corners.

Moment maps of molecular gas can provide insight into the morphology and
kinematics of the gas. Moment maps were generated using the 3D-Based Analysis
of Rotating Objects via Line Observations (3D-Barolo; Di Teodoro & Fraternali,
2015), which employs the Duchamp algorithm to identify pixels above a primary
signal-to-noise ratio (S/N) and agglomerates neighboring pixels above a
secondary S/N. The initial moment maps were created by running 3D-Barolo on
the non-primary beam (PB) attenuation corrected image cubes, and the final
moment maps were derived by multiplying the PB attenuation corrected image
cubes by the mask image cubes created during the generation of the initial
moment maps. For this analysis, the primary and secondary S/N were set at 5
and 3, respectively.

In Figure 4, the moment maps of the CO (1–0), CO (3–2), HCN (1–0), HCO+ (1–0),
and HCO+ (3–2) lines are shown from top to bottom. The total intensity (moment
0), velocity field (moment 1), and velocity dispersion (moment 2) maps are
displayed in the left, middle, and right panels, respectively. The black
contours indicate the 99 GHz continuum, with the BH located at the local
maximum between two major blobs of star formation. The moment maps of the
molecular gas lines indicate that there is a significant amount of molecular
gas beyond the continuum regions traced by the black contours. These extended
features are consistent with previous observations (Baas et al., 1994; Beck et
al., 2018; Imara & Faesi, 2019). Total intensity maps of the molecular gas
(left panels of Figure 4) show two prominent regions corresponding to the
strong continuum blobs of Figure 2. The maximum intensities of the molecular
gas in the left panels of Figure 4 coincide with the maximum positions in the
continuum. The BH is located between the two prominent starburst regions in
the molecular gas distributions, forming a spur-like structure visible in the
high-resolution CO data. The detection of HCO+ (1–0), HCO+ (3–2), and HCN
(1–0) above S/N $>$ 5 indicate the existence of significant dense gas (with
critical density $\gtrsim~{}10^{6}$ cm-3) in the vicinity of the BH.

####  3.2.2 Position-Velocity Diagram

Figure 5: PVDs of CO J$=$3–2 (panel A) and CO J$=$1–0 (panel B) are presented
with the weighted mean velocity in each position bin above 3$\sigma$ (red
squares) and the velocity gradient (black lines). The light blue contours
indicate 3, 10, and 50$\sigma$ in CO (3–2) and 3, 5, and 10$\sigma$ in CO
(1–0), respectively. The dotted vertical lines present the regions estimating
the velocity gradient, e.g., -0$\farcs$4 $<x<$ 0$\farcs$8 for CO (3–2) and
-0$\farcs$5 $<x<$ 0$\farcs$85 for CO (1–0).

The position-velocity diagram (PVD) serves as a valuable tool for elucidating
the kinematics of gas. The ionized gas in the central region of He 2–10
exhibits a sinusoidal motion in the PVD, which is well fit by the bipolar
outflow model presented in Schutte & Reines (2022). Schutte & Reines (2022)
used HST/STIS spectroscopy, employing a long slit centered on the BH with a
width of 0$\farcs$2 and a position angle of 100∘ aligned with the ionized
filament seen in H$\alpha$ imaging. In our investigation, we look for a
similar pattern in the velocities of the molecular gas by generating PVDs
along a line with a position angle of 100∘, centered on the BH, from our CO
(3–2) and CO (1–0) image cubes. We extracted PVDs along the line in bins
corresponding to the synthesized beam sizes of CO (3–2) and CO (1–0) image
cubes, i.e., 0$\farcs$3 and 0$\farcs$45, respectively. However, these values
are more than 3$\times$ larger than the spatial resolution of the HST/STIS
observations ($\sim 0\farcs 1$).

Figure 5 illustrates the resulting PVDs for CO (3–2) in panel (A) and CO (1–0)
in panel (B). Unlike the ionized gas observed by Schutte & Reines (2022), no
sinusoidal motions were discerned in either PVD. This may be attributed to the
coarser angular resolution of our observations, leading to a beam smearing
effect that obscures outflow signals amid adjacent background emissions.
Alternatively, the absence of sinusoidal motion in the molecular gas could
signify a genuine lack of such motion, suggesting a potential insensitivity of
the molecular gas to the BH outflow, particularly if the BH outflow lacks
sufficient strength (Koudmani et al., 2022).

We calculated the velocity gradient from the PVD of the CO (3–2) line in panel
(A) of Figure 5, leveraging the smaller beam size of the CO (3–2) image cube
compared to that of CO (1–0). We restrict this analysis to the gas in the
central region near the BH between -0$\farcs$4 ($-18$ pc) to 0$\farcs$8 ($35$
pc) with velocities in the range $\sim$870–900 km s-1. The resulting weighted
mean velocity (red squares in Figure 5), determined above 3$\sigma$ with a bin
width of 0$\farcs$3 (where $\sigma=1.0$ mJy beam-1 (1.693 km s-1)-1, yielded a
velocity gradient of $-5.8$ km s-1 ″-1. We repeated this analysis for CO
(1–0), constructing a PVD along a line with a position angle of 100∘ and a
width of 0$\farcs$45 (matching the beam size of CO (1–0)), resulting in a
similar velocity gradient of $-4.6$ km s-1 ″-1).

Our derived velocity gradient is shallower than the $-7.8$ km s-1 ″-1 reported
by Imara & Faesi (2019). The discrepancy can be attributed to differences in
i) the length of the position axis used to estimate the velocity gradient, ii)
the position angle, and iii) the bin size used to trace gas velocity. Imara &
Faesi (2019) estimated the velocity gradient over a range of -1$\farcs$7
($\sim$70 pc) to 1$\farcs$7 along a line with a position angle of 143∘ and a
width of 4″. In our PVD (Figure 5), molecular gas at approximately 1$\farcs$5
is associated with the western star-forming region, contributing to a steeper
velocity gradient due to its high velocity. The differences in position angle
and bin size introduce additional inconsistencies, particularly given the
asymmetrical distribution of molecular gas around the BH.

####  3.2.3 Molecular Gas in the Vicinity of the BH

Figure 6: Moment maps of the molecular gas species observed at high angular
resolution including CO (1–0), CO (3–2), and HCO+ (3–2) (top to bottom), where
the moment maps are generated with a restrictive velocity range of 852.1–904.6
km s-1 (see §3.2). The total intensity, velocity field, and velocity
dispersion are presented in the left, middle, and right panels, respectively.
The small black contour indicates the position of BH traced by the high
brightness temperature radio core detected with the LBA (Reines & Deller,
2012). The outer black ellipse defines the region we refer to as the BH
vicinity, which has a near uniform velocity that is consistent with the
systemic velocity of the galaxy, $V_{sys}=873$ km s-1 (Kobulnicky et al.,
1995).

In this work, we are primarily interested in examining the molecular gas in
the vicinity of the BH333The ALMA observations do not have sufficient angular
resolution to probe the BH sphere of influence. Given the BH mass of $\rm\sim
10^{6}\;M_{\odot}$ (Reines et al., 2016; Riffel, 2020), the sphere of
influence is $\sim 4$ pc (0$\farcs$1), which is at least three times smaller
than the angular resolution of our observations.. We define this region (i.e.,
the BH vicinity) as follows. There is a distinct velocity gas blob near the BH
in HCO+ (3–2), seen in the bottom middle panel of Figure 4. The velocity range
of this blob is 852.1–904.6 km s-1 in HCO+(3–2). To help isolate the gas in
the vicinity of the BH, both spatially and in velocity, we generated moment
maps of the molecular gas tracers with higher angular resolution (CO (1–0), CO
(3–2), and HCO+ (3–2)) within this restricted velocity range (Figure 6).

The resulting moment maps reveal an elongated CO gas structure with a distinct
velocity from the surrounding regions, as seen in the middle panels of Figure
6. This distinct velocity region is defined as the BH vicinity, represented by
the black ellipse in Figure 6. The ellipse is centered at RA 08h 36m 15.108s
and Dec $-26$∘ 24′ 34.343″, with a size of 2$\farcs$97$\times$0$\farcs$79
($\approx$ 128 pc$\times$34 pc) and a position angle of 80.34∘ as presented in
Table 3. The moment maps in Figure 6 clearly demonstrate that the BH (small
black contour) is located north of the BH vicinity (black ellipse) in
molecular gases such as CO and HCO+ (3–2). In other words, the BH is not at
the center of the molecular gas structure. The lack of the molecular gas to
the north of the BH is likely a result of star formation feedback occurring in
the young massive star clusters with ages of $\sim 4-5$ Myrs (Beck et al.,
2018), as seen in the H${\alpha}$ image in Figure 7. This could help explain
why the BH is accreting at such a low Eddington ratio ($\sim 10^{-6}$; Reines
et al. 2016). Furthermore, the asymmetrical distribution of molecular gas
around the massive BH resembles the central molecular zone of the Milky Way
(Bally et al., 1988) and the central region of NGC 253 (Meier et al., 2015;
Krieger et al., 2019; Levy et al., 2022). In the cases of the Milky Way and
NGC 253, the origin of this asymmetry is attributed to feedback from star
formation in their centers (Henshaw et al., 2022; Krieger et al., 2019).

We have measured the velocity ranges within the BH vicinity using three
different molecular lines: CO J=1–0, CO J=3–2, and HCO+ J=3–2. The velocity
ranges were found to be 865–899 km s-1, 871–904 km s-1, and 871–904 km s-1,
respectively. Excluding the borders, the velocity ranges for CO J=1–0 and
J=3–2 were reduced to 871–886 km s-1 and 871–885 km s-1, respectively. The
mean velocities within the BH vicinity are 880$\pm$5 (CO J=1–0), 880$\pm$6 (CO
J=3–2), and 875$\pm$6 km s-1 (HCO+ J=3–2), and are therefore consistent with
one another within the uncertainties. For reference, the systemic velocity of
the galaxy is 873 km s-1 (Kobulnicky et al., 1995). The velocity field maps do
not show any significant rotation in the region defined as the BH vicinity.
The velocity dispersion map of CO J=3–2 shows some high velocity dispersion
along the SE and NW directions with $\Delta\sigma\approx 7$ km s-1. The high
velocity dispersion ($\Delta\sigma\approx 50$ km s-1) of stellar motions near
the BH observed by Riffel (2020) is not seen in the molecular gas.

Table 3: Properties of molecular gas in the vicinity of the BH (black ellipse in Figure 6) Parameter | Value
---|---
Center coordinate | RA 08h 36m 15.108s, Dec $-26$∘ 24′ 34.343″
Major & minor axes | 2$\farcs$97$\times$0$\farcs$79 ($\approx$128 pc$\times$34 pc)
Position angle | 83.34∘
Mean velocity | 880.2 $\pm$ 4.7 km s-1 in CO (1–0)
880.3 $\pm$ 5.9 km s-1 in CO (3–2)
875.3 $\pm$ 6.2 km s-1 in HCO+ (3–2)
H2 mass | (0.21–1.53)$\times 10^{6}$ M⊙

The molecular gas mass within the BH vicinity is estimated using the scaling
relation between molecular gas mass and CO luminosity, $\rm
M_{H_{2}}=\alpha_{CO}L^{{}^{\prime}}_{CO}$. The conversion factor
$\alpha_{CO}$ depends on the metallicity (Bolatto et al., 2013, and references
therein). Cresci et al. (2017) conducted spectroscopic observations with MUSE
to determine the metallicity of He 2–10, with a spatial sampling size of
0$\farcs$2$\times$0$\farcs$2\. They created a map of metallicity distribution,
which revealed significant variation across the galaxy, largely due to change
in ionization parameters. Through simultaneous line-fitting on multiple lines,
they found that the metallicity in the central region was super-solar, with
$\rm 12+[O/H]\approx 9.0$, corresponding to $\rm Z=0.044$. Based on the
relationship between $\alpha_{\rm CO}$ and metallicity outlined in Bolatto et
al. (2013), we choose a range of $\alpha_{\rm CO}$ values between 0.6 and 4.3
$M_{\odot}$ (K km s-1 pc2)-1, where $\rm\alpha_{\rm CO}=4.3$ is the value for
the Milky Way. Our range of $\alpha_{\rm CO}$ also includes the widely
accepted value for the centers of starburst galaxies with solar metallicities,
$\alpha_{\rm CO}\approx 1.1$ M⊙ (K km s-1 pc2)-1 (Leroy et al., 2015; Krieger
et al., 2019). These considerations ensure that our analysis is consistent
with established standards in the field. Our estimated CO luminosity is $\rm
L^{{}^{\prime}}_{CO}=3.56\times 10^{5}\;K\;km\;s^{-1}\;pc^{2}$, which when
combined with the $\rm\alpha_{\rm CO}$ range, results in an estimated
molecular gas mass of $\rm M_{H_{2}}=(0.21-1.53)\times 10^{6}$ M⊙.

####  3.2.4 CO Line Ratios

Figure 7: Left: Pixel-by-pixel line ratio map between CO (3–2) and CO (1–0),
$R_{31}$. Right: HST H${\alpha}$ plus continuum image with the same field of
view. The BH position is marked with black contours from the LBA 1.4 GHz
observations by Reines & Deller (2012), and the region defined as the BH
vicinity is indicated with the orange dashed ellipse. Regions 5 and 6 from
Schutte & Reines (2022) are presented as rectangles with the HST slit position
for the ionized gas outflow (Schutte & Reines, 2022). The highest line ratios
are observed at the interface between the BH outflow and the eastern star
forming region, as well as region 5 to the west of the BH. This implies the
cold molecular gas in these regions is shocked by the bipolar outflow,
supporting the case for positive BH feedback.

Line ratios between different molecular gas tracers are good indicators of gas
excitation (Privon et al., 2015). The HCN (1–0) and HCO+ (1–0) line image
cubes have large beam sizes, and the HCO+ (3–2) line image cube has a
relatively low signal-to-noise ratio, so we restricted our analysis to the
line ratio between CO J=1–0 and J=3–2, which have higher angular resolutions
and better sensitivities. Given that the CO (3–2) line emission exhibits a
more extended distribution compared to the CO (1–0) line emission, the
distribution of the CO line ratio is constrained by the CO (1–0) line
distribution. The line ratio is defined as $\rm
R_{31}=L_{CO(3-2)}/L_{CO(1-0)}=(\nu_{32}/\nu_{10})^{-2}(I_{CO(3-2)}/I_{CO(1-0)})$,
where $\rm I_{CO}$ is the velocity integrated line intensity in units of Jy km
s-1. Both image cubes were convolved to have the same synthesized beam size of
0$\farcs$45$\times$0$\farcs$30 using the CASA task imsmooth. We present the
pixel-by-pixel CO line ratio map in the left panel and an HST H$\alpha$ image
(e.g., Schutte & Reines, 2022) in the right panel of Figure 7. Our line ratio
map of the central region of the galaxy indicates highly excited CO molecular
gas as evidenced by $\rm R_{31}\gtrsim 1$, exceeding the typical range
observed in the Milky Way, $\rm 0.3<R_{31}<0.5$ (Leroy et al., 2009; Carilli &
Walter, 2013). This result aligns with the findings of Beck et al. (2018), who
reported a high flux ratio of $\rm S_{CO(3-2)}/S_{CO(1-0)}=6.1$ for this
galaxy.

####  3.2.5 Shocked Molecular Clouds from the BH Outflow

The discovery of a low-velocity bipolar outflow originating from the BH was
reported in ionized gas (Schutte & Reines, 2022). The precession angle and
frequency of this outflow were estimated to be in the range of $2.4-6.1$∘ and
$3.0-7.5$ revolutions per Myr, respectively. Notably, Schutte & Reines (2022)
identified several impacts of the BH outflow, including the stimulation of
star formation in the eastern star-forming region, as well as in the dark
cloud (region 5) and the dense clouds within the western star-forming region
(region 6). The left panel of Figure 7 reveals that the highest $\rm R_{31}$
line ratios are observed at the interface between the known ionized gas
outflow from the BH and the eastern star forming region (Schutte & Reines,
2022), as depicted in the H$\alpha$ image (right panel of Figure 7). The line
ratio at this interface is $\rm R_{31}=2.7\pm 0.2$, while the region defined
as the BH vicinity excluding this interface has $R_{31}=1.3\pm 0.1$.
Additionally, we detect a high line ratio ($\rm R_{31}=2.5\pm 0.7$) to the
west of the BH at the boundary of a dark cloud of gas and dust (region 5 in
Figure 7). This region was previously identified by Schutte & Reines (2022)
and exhibits double-peaked optical emission lines, likely due to the other
side of the bipolar BH outflow intercepting dense clouds and pushing them
laterally (as is the case with the eastern interface region). The CO line
ratio at region 5 is comparable to that of interface between the outflow and
the eastern star forming region, consistent with this scenario. These findings
suggest that the molecular gas impacted by the BH outflow to the east and west
have both high density and high kinetic temperature (Carilli & Walter, 2013).
It is noteworthy that the elevated $\rm R_{31}$ may be an edge effect limited
by low flux densities of CO (1–0). Nevertheless, certain regions at the
peripheries, particularly in the northeast and south, exhibit low $\rm R_{31}$
values. The presence of these regions implies that the the pixel-by-pixel $\rm
R_{31}$ values are influenced by the relative distributions of CO (1–0) and CO
(3–2) due to some physical origins, rather than being solely dictated by low
CO (1–0) flux densities. Moreover, our derived $\rm R_{31}$ values are robust,
exhibiting similar sensitivities, e.g., $\sigma=$1.03 mJy beam-1 (1.693 km
s-1)-1 for CO (1–0) and 1.00 mJy beam-1 (1.693 km s-1)-1 for CO (3–2),
respectively. Therefore, the elevated $\rm R_{31}$ values observed in this
study are not solely constrained by low CO (1–0) flux densities.

Multiple studies have investigated the occurrence of elevated CO line ratios
in different galaxies. García-Burillo et al. (2014) and Viti et al. (2014)
examined the CO line ratio in NGC 1068 and demonstrated that the ratio between
J=3–2 and J=1–0 vary significantly between the circumnuclear disk around the
AGN ($R_{31}=2.7$) and the starburst ring ($R_{31}=1.2$). Furthermore, they
observed higher CO line ratios in more active star forming regions compared to
less active ones ($R_{31}=0.7-1.0$). Their analysis suggested that the
elevated CO line ratio in the circumnuclear disk is a result of the AGN-driven
outflow. The CO line ratio in IC 5063 was extensively studied by Dasyra et al.
(2016) and Oosterloo et al. (2017). They identified excited molecular gas
along the radio jet and noted a higher CO line ratio in the central region
compared to the outer region. For instance, Dasyra et al. (2016) conducted a
calculation of the line ratio between CO (4–3) and CO (2–1), revealing a value
of 12 in the robust jet-cloud interacting regions compared to an average of 5.
Oosterloo et al. (2017) demonstrated that $R_{31}>1$ along the jets as opposed
to 0.17 in the disk. Their study further indicated that the most elevated CO
line ratio is proximate to the two radio lobes, implicating the shock
generated by the outflow as the causative factor (Morganti et al., 2021). More
recently, Audibert et al. (2023) identified a high CO line ratio between
J$=3-2$ and J$=2-1$, specifically with $T_{CO(3-2)}/T_{CO(2-1)}=0.8$ in the
region perpendicular to the radio jet, in contrast to the value of 0.4 in
other regions within the type-2 QSO Teacup galaxy. Their findings underscored
that this augmented CO line ratio is a consequence of the jet-driven shock
interacting with the molecular gas.

Our observations align with the results reported by García-Burillo et al.
(2014), where the $R_{31}$ values stand at 2.7 in the circumnuclear disk, 1.2
in starburst regions, and 0.4$-$0.5 in the halo. Consequently, the regions
exhibiting high $R_{31}$ values in our study can be attributed to shock
heating induced by the BH outflow.

####  3.2.6 The Western Star Forming Region

While prominent at radio wavelengths (e.g., Figure 2), the star-forming region
to the west of the BH appears to suffer from significant extinction at optical
wavelengths (e.g., Figures 8 and 7). Schutte & Reines (2022) estimate a
stellar age of $<$3 Myr for region 6 (e.g., see Figure 7) based on the
equivalent width of the H$\alpha$ emission line. This is consistent with a
scenario in which the star clusters have not had enough time to clear away
their birth material, leading to high levels of extinction. We quantify their
assertion with the molecular gas distribution.

Figure 8: We present the channel maps of the CO (3–2) line within the velocity
range of 864.1–887.8 km s-1. In order to highlight the faint emission beyond
the BH vicinity (indicated by the dashed orange line), we have masked the
majority of CO gas. The 99 GHz continuum emission is depicted by black
contours at levels of 50, 100, 200, and 300 $\mu$Jy beam-1. The BH is
represented by a solid black ellipse, corresponding to the central maximum of
the 1.4 GHz continuum observed in LBA (Reines & Deller, 2012). The diagonal
dotted lines illustrate the spectroscopic slit used in the HST observations,
with a width of 0$\farcs$2\. Additionally, the regions 5 (left) and 6 (right)
previously discussed by Schutte & Reines (2022) are depicted by red boxes. The
ALMA beam size is provided in the bottom-right corner of each panel. The
bottom row exhibits the total intensity map of the CO (3–2) line on the left,
followed by the HST I-band image in the middle, and the HST
H$\alpha$+continuum image on the right. In the H$\alpha$+continuum image
(panel r), white contours represent the 99 GHz continuum map, cyan dotted
lines depict the slit, and green boxes highlight regions 5 and 6 for improved
visibility.

Figure 8 presents channel maps of the CO (3–2) line within velocities of
864–888 km s-1 in panels (a) to (o), a total intensity map in panel (p), an
HST I-band image in panel (q), and an HST H$\alpha$ plus continuum image in
panel (r). The total intensity map was produced using 3D-Barolo (Di Teodoro &
Fraternali, 2015) with a primary S/N of 5 and a second S/N of 3 on the image
cube shown in panels (a)–(o). In the channel maps and total intensity map, the
bulk of the gas was masked to enhance contrast and highlight fainter emission.
The red boxes (white boxes in panel r) in each panel of Figure 8 indicate
regions 5 and 6 defined by Schutte & Reines (2022), respectively. Figure 8
demonstrates that there is a significant amount of molecular gas between
regions 5 and 6, implying that the molecular gas is abundant enough to absorb
the emission from young stellar populations embedded in the gas.

We estimate the optical extinction in this region as follows. The total
intensity of the CO (1–0) line within the span of regions 5 and 6 is measured
as $\rm W_{CO(1-0)}=788.2\pm 57.4$ K km s-1, resulting in an estimated
molecular gas column density of $\rm N(H_{2})=(2.21\pm 0.16)\times
10^{22}-(1.58\pm 0.12)\times 10^{23}$ cm-2. This estimation is derived from
the mean velocity-integrated flux density ratio of CO (3–2) to CO (1–0),
represented as $\rm I_{CO(3-2)}/I_{CO(1-0)}=12.54$ in this region. The
conversion factor $\rm X_{CO}=(0.28-2)\times 10^{20}$ cm-2 (K km s-1)-1,
corresponding to $\rm\alpha_{CO}=0.6-4.3$ $M_{\odot}$ (K km s-1 pc2)-1, has
been employed for this calculation. Kobulnicky et al. (1995) conducted
observations of atomic hydrogen (H I) using the VLA at D-, C-, and B-arrays.
In this region, the total intensity of H I was 1.52$\pm$0.05 Jy beam-1 km s-1,
leading to a corresponding column density N(H I) at
(1.51$\pm$0.05)$\times$1018 cm-2, considering a synthesized beam size of
31″$\times$30″. Assuming a uniform distribution of H I, the optical extinction
can be computed using the relationship between N(H) and AV (Zhu et al., 2017),
$\rm A_{V}=4.81\times 10^{-22}(N(HI)+2N(H_{2}))$. This calculation yields a
range of $\rm A_{V}\sim 21-152$ mag. If the H I is more concentrated in this
region, these values would represent lower limits. On the other hands, the
optical extinction is estimated to be $\rm A_{V}\sim 7-51$ for the eastern
star forming region, supporting that the optically dark region is
predominantly observed between region 5 and 6.

The 251 and 340 GHz images presented in Figure 2 showed that the compact
features of region 6 align with peaks in CO (3–2) and HCO+ (3–2) shown in
Figure 4. This observation is in line with the characteristics of bright and
compact embedded clusters seen in other young clusters, e.g., NGC 253 (Ando et
al., 2017; Leroy et al., 2018; Levy et al., 2021; Mills et al., 2021) and NGC
4945 (Emig et al., 2020). To quantify our findings, we focused on the radio
continuum blob located in the western star-forming region, corresponding to
region 6 as defined in Schutte & Reines (2022). Dust and gas mass estimations
were carried out within this region. For the calculation of dust mass, we
utilized the 340 GHz flux density and employed the equation:

$M_{dust}={{S_{\nu}D^{2}}\over{\kappa_{\nu}B_{\nu}(T_{d})}}$ (1)

(Launhardt & Henning, 1997). We adopted a dust absorption coefficient of
$\kappa_{\nu}=1.3$ cm2 g-1 at 340 GHz (Costa et al., 2021) and the dust
temperature of 29 K (Santangelo et al., 2009). The measured flux density at
340 GHz is $S_{340GHz}=610.6\pm 30.0$ $\mu$Jy, resulting in a calculated dust
mass of $(2.37\pm 0.12)\times 10^{3}$ M⊙. The associated molecular mass is
$M_{\rm H_{2}}=(7.79\pm 0.15)\times 10^{4}-(5.65\pm 0.11)\times 10^{5}$ Modot
depending on the $\alpha_{CO}$. Consequently, the dust to molecular gas mass
ratio falls within the range of 1:33–238, which is consistent with the typical
value of Milky Way and other starburst galaxies, 1:120 (Wilson et al., 2008).

The existence of abundant molecular gas and high optical extinction value in
the optically dark regions (regions 5 and 6) indicate that the molecular gas
in these regions is compressed. In He 2–10, the shock generated by the BH
outflow emerges as the most plausible source for the compression of molecular
gas, contributing to the formation of embedded star clusters. This is an
example of the “positive feedback” of the BH outflow.

The asymmetric distributions of molecular gas and embedded star formation
provides a key to understanding the evolution of He 2–10. Given that He 2–10
has experienced a tidal interaction (Kobulnicky et al., 1995), the asymmetric
distributions of molecular gas result from the tidal interaction by disturbing
the galaxy potential. This asymmetric distribution of molecular gas might
introduce the asymmetric embedded star formations via the asymmetric gas
inflow, similar to the Central Molecular Zone of the Milky Way (Battersby et
al., 2020; Henshaw et al., 2022).

##  4 Conclusions

We re-imaged the ALMA archival data of He 2–10 using newly developed
algorithms and studied properties of the continuum and molecular gas,
including CO (1–0), CO (3–2), HCN (1–0), HCO+ (1–0), and HCO+ (3–2), focusing
on the massive BH and its vicinity. We re-analyzed the radio–sub-mm SED of the
BH with updated flux densities and investigated the kinematics and morphology
of molecular gas around the BH. Our main results are summarized below.

  * •

The radio–sub-mm SED of the BH from $1.4-340$ GHz is dominated by synchrotron
emission with a best-fit model of $S_{\nu}=(1349\pm 217)\nu^{-0.52\pm 0.06}$
$\mu$Jy (Figure 3).

  * •

We find an elongated molecular gas structure ($\approx 128~{}{\rm pc}\times
34~{}{\rm pc}$) in the vicinity of the BH with a distinct velocity component
in CO (1–0), CO (3–2), and HCO+ (3–2) from the surrounding regions (Figure 6).
The mass of the molecular gas structure in the vicinity of the BH is $\rm
M_{H_{2}}=(0.21-1.53)\times 10^{6}\;M_{\odot}$ for a range of metallicity-
dependent values of $\alpha_{CO}$. While the position of the BH is within this
region, it is significantly offset from the peak intensity, which could
explain why the BH is accreting at a low rate.

  * •

The ratio of CO (3–2) to CO (1–0) shows that the CO gas is highly excited in
the BH vicinity (Figure 7). The highest excitation values are located at the
interface between the BH outflow and eastern region of triggered star
formation identified by Schutte & Reines (2022). Our results suggest that the
molecular gas is being compressed by the shock produced by the outflow from
the BH, providing additional evidence for positive BH feedback.

  * •

The high excitation of the CO line ratio is also detected to the west of the
BH, suggesting the other side of the bipolar flow is shocking dense molecular
clouds in the western star forming region as well. We find high molecular gas
column densities (${\rm N(H_{2})\sim 10^{21-22}}~{}{\rm cm}^{-2}$) in this
region, as well as high levels of extinction ($A_{V}\gtrsim 2-11$ mag). This
is consistent with a scenario in which the infant star clusters with stellar
ages $\lesssim$ 3 Myr (Schutte & Reines, 2022) have not yet had enough time to
destroy or disperse the dense molecular clouds from which they formed.

The study presented in this paper constitutes the first investigation of the
molecular gas properties in the immediate vicinity of the massive BH in He
2–10. While our results provide additional support for positive BH feedback,
future observations with higher angular resolution are needed to resolve the
sphere of influence of the BH and trace any potential molecular gas outflow
directly. Furthermore, analyzing the higher transitions of the CO molecular
gas and applying models for the photodissociation region (PDR) and/or X-ray
dominated region (XDR) will provide additional information on physical
parameters such as density, temperature, and outflow power, ultimately
yielding more insight into the impact of BH feedback on nearby molecular gas.

## Acknowledgement

The authors express their gratitude to the anonymous reviewer for their
valuable contributions in enhancing the quality of this paper. We are grateful
to John Hibbard, Loreto Barcos-Munoz, Mark Lacy, Allison Costa, George Privon
at NRAO, Seth Kimbrell, and Sheyda Salehirad at Montana State University for
valuable discussions. A.E.R. acknowledges support for this work provided by
NASA through EPSCoR grant number 80NSSC20M0231. This paper makes use of the
following ALMA data: ADS / JAO.ALMA #2011.0.00348.S, ADS / JAO.ALMA
#2012.1.00413.S, ADS / JAO.ALMA #2015.1.01569.S, ADS / JAO.ALMA
#2016.1.00027.S, ADS / JAO.ALMA #2016.1.00492.S, ADS / JAO.ALMA
#2019.1.01641.S. ALMA is a partnership of ESO (representing its member
states), NSF (USA) and NINS (Japan), together with NRC (Canada), MOST and
ASIAA (Taiwan), and KASI (Republic of Korea), in cooperation with the Republic
of Chile. The Joint ALMA Observatory is operated by ESO, AUI/NRAO and NAOJ.
This research has made use of NASA’s Astrophysics Data System Bibliographic
Services.

## References

  * Allen et al. (1976) Allen, D. A., Wright, A. E., & Goss, W. M. 1976, MNRAS, 177, 91, doi: 10.1093/mnras/177.1.91
  * Ando et al. (2017) Ando, R., Nakanishi, K., Kohno, K., et al. 2017, ApJ, 849, 81, doi: 10.3847/1538-4357/aa8fd4
  * Audibert et al. (2023) Audibert, A., Ramos Almeida, C., García-Burillo, S., et al. 2023, A&A, 671, L12, doi: 10.1051/0004-6361/202345964
  * Baas et al. (1994) Baas, F., Israel, F. P., & Koornneef, J. 1994, A&A, 284, 403
  * Bally et al. (1988) Bally, J., Stark, A. A., Wilson, R. W., & Henkel, C. 1988, ApJ, 324, 223, doi: 10.1086/165891
  * Battersby et al. (2020) Battersby, C., Keto, E., Walker, D., et al. 2020, ApJS, 249, 35, doi: 10.3847/1538-4365/aba18e
  * Bayet et al. (2004) Bayet, E., Gerin, M., Phillips, T. G., & Contursi, A. 2004, A&A, 427, 45, doi: 10.1051/0004-6361:20035614
  * Bean et al. (2022) Bean, B., Bhatnagar, S., Castro, S., et al. 2022, PASP, 134, 114501, doi: 10.1088/1538-3873/ac9642
  * Beck et al. (2018) Beck, S. C., Turner, J. L., & Consiglio, S. M. 2018, ApJ, 867, 165, doi: 10.3847/1538-4357/aae481
  * Bolatto et al. (2013) Bolatto, A. D., Wolfire, M., & Leroy, A. K. 2013, ARA&A, 51, 207, doi: 10.1146/annurev-astro-082812-140944
  * Carilli & Walter (2013) Carilli, C. L., & Walter, F. 2013, ARA&A, 51, 105, doi: 10.1146/annurev-astro-082812-140953
  * Comrie et al. (2021) Comrie, A., Wang, K.-S., Hsu, S.-C., et al. 2021, CARTA: The Cube Analysis and Rendering Tool for Astronomy, 2.0.0, Zenodo, Zenodo, doi: 10.5281/zenodo.4905459
  * Costa et al. (2021) Costa, A. H., Johnson, K. E., Indebetouw, R., et al. 2021, ApJ, 918, 76, doi: 10.3847/1538-4357/ac0e93
  * Cresci et al. (2017) Cresci, G., Vanzi, L., Telles, E., et al. 2017, A&A, 604, A101, doi: 10.1051/0004-6361/201730876
  * Dasyra et al. (2016) Dasyra, K. M., Combes, F., Oosterloo, T., et al. 2016, A&A, 595, L7, doi: 10.1051/0004-6361/201629689
  * Di Teodoro & Fraternali (2015) Di Teodoro, E. M., & Fraternali, F. 2015, MNRAS, 451, 3021, doi: 10.1093/mnras/stv1213
  * Emig et al. (2020) Emig, K. L., Bolatto, A. D., Leroy, A. K., et al. 2020, ApJ, 903, 50, doi: 10.3847/1538-4357/abb67d
  * García-Burillo et al. (2014) García-Burillo, S., Combes, F., Usero, A., et al. 2014, A&A, 567, A125, doi: 10.1051/0004-6361/201423843
  * Greene et al. (2020) Greene, J. E., Strader, J., & Ho, L. C. 2020, ARA&A, 58, 257, doi: 10.1146/annurev-astro-032620-021835
  * Greisen (2003) Greisen, E. W. 2003, in Astrophysics and Space Science Library, Vol. 285, Information Handling in Astronomy - Historical Vistas, ed. A. Heck, 109, doi: 10.1007/0-306-48080-8_7
  * Henshaw et al. (2022) Henshaw, J. D., Krumholz, M. R., Butterfield, N. O., et al. 2022, MNRAS, 509, 4758, doi: 10.1093/mnras/stab3039
  * Imanishi et al. (2007) Imanishi, M., Nakanishi, K., Tamura, Y., Oi, N., & Kohno, K. 2007, AJ, 134, 2366, doi: 10.1086/523598
  * Imara & Faesi (2019) Imara, N., & Faesi, C. M. 2019, ApJ, 876, 141, doi: 10.3847/1538-4357/ab16cc
  * Johnson et al. (2018) Johnson, K. E., Brogan, C. L., Indebetouw, R., et al. 2018, ApJ, 853, 125, doi: 10.3847/1538-4357/aa9ff8
  * Johnson & Kobulnicky (2003) Johnson, K. E., & Kobulnicky, H. A. 2003, ApJ, 597, 923, doi: 10.1086/378585
  * Kobulnicky et al. (1995) Kobulnicky, H. A., Dickey, J. M., Sargent, A. I., Hogg, D. E., & Conti, P. S. 1995, AJ, 110, 116, doi: 10.1086/117500
  * Kormendy & Ho (2013) Kormendy, J., & Ho, L. C. 2013, ARA&A, 51, 511, doi: 10.1146/annurev-astro-082708-101811
  * Koudmani et al. (2022) Koudmani, S., Sijacki, D., & Smith, M. C. 2022, MNRAS, 516, 2112, doi: 10.1093/mnras/stac2252
  * Krieger et al. (2019) Krieger, N., Bolatto, A. D., Walter, F., et al. 2019, ApJ, 881, 43, doi: 10.3847/1538-4357/ab2d9c
  * Launhardt & Henning (1997) Launhardt, R., & Henning, T. 1997, A&A, 326, 329
  * Leroy et al. (2009) Leroy, A. K., Walter, F., Bigiel, F., et al. 2009, AJ, 137, 4670, doi: 10.1088/0004-6256/137/6/4670
  * Leroy et al. (2015) Leroy, A. K., Bolatto, A. D., Ostriker, E. C., et al. 2015, ApJ, 801, 25, doi: 10.1088/0004-637X/801/1/25
  * Leroy et al. (2018) —. 2018, ApJ, 869, 126, doi: 10.3847/1538-4357/aaecd1
  * Levy et al. (2021) Levy, R. C., Bolatto, A. D., Leroy, A. K., et al. 2021, ApJ, 912, 4, doi: 10.3847/1538-4357/abec84
  * Levy et al. (2022) —. 2022, ApJ, 935, 19, doi: 10.3847/1538-4357/ac7b7a
  * Meier et al. (2001) Meier, D. S., Turner, J. L., Crosthwaite, L. P., & Beck, S. C. 2001, AJ, 121, 740, doi: 10.1086/318782
  * Meier et al. (2015) Meier, D. S., Walter, F., Bolatto, A. D., et al. 2015, ApJ, 801, 63, doi: 10.1088/0004-637X/801/1/63
  * Mills et al. (2021) Mills, E. A. C., Gorski, M., Emig, K. L., et al. 2021, ApJ, 919, 105, doi: 10.3847/1538-4357/ac0fe8
  * Mohan & Rafferty (2015) Mohan, N., & Rafferty, D. 2015, PyBDSF: Python Blob Detection and Source Finder, Astrophysics Source Code Library, record ascl:1502.007.  http://ascl.net/1502.007
  * Morganti et al. (2021) Morganti, R., Oosterloo, T., Murthy, S., & Tadhunter, C. 2021, Astronomische Nachrichten, 342, 1135, doi: 10.1002/asna.20210037
  * Nguyen et al. (2014) Nguyen, D. D., Seth, A. C., Reines, A. E., et al. 2014, ApJ, 794, 34, doi: 10.1088/0004-637X/794/1/34
  * Oosterloo et al. (2017) Oosterloo, T., Raymond Oonk, J. B., Morganti, R., et al. 2017, A&A, 608, A38, doi: 10.1051/0004-6361/201731781
  * Privon et al. (2015) Privon, G. C., Herrero-Illana, R., Evans, A. S., et al. 2015, ApJ, 814, 39, doi: 10.1088/0004-637X/814/1/39
  * R Core Team (2022) R Core Team. 2022, R: A Language and Environment for Statistical Computing, R Foundation for Statistical Computing, Vienna, Austria.  https://www.R-project.org/
  * Rau & Cornwell (2011) Rau, U., & Cornwell, T. J. 2011, A&A, 532, A71, doi: 10.1051/0004-6361/201117104
  * Reines (2022) Reines, A. E. 2022, Nature Astronomy, 6, 26, doi: 10.1038/s41550-021-01556-0
  * Reines & Deller (2012) Reines, A. E., & Deller, A. T. 2012, ApJ, 750, L24, doi: 10.1088/2041-8205/750/1/L24
  * Reines et al. (2016) Reines, A. E., Reynolds, M. T., Miller, J. M., et al. 2016, ApJ, 830, L35, doi: 10.3847/2041-8205/830/2/L35
  * Reines et al. (2011) Reines, A. E., Sivakoff, G. R., Johnson, K. E., & Brogan, C. L. 2011, Nature, 470, 66, doi: 10.1038/nature09724
  * Reines & Volonteri (2015) Reines, A. E., & Volonteri, M. 2015, ApJ, 813, 82, doi: 10.1088/0004-637X/813/2/82
  * Riffel (2020) Riffel, R. A. 2020, MNRAS, 494, 2004, doi: 10.1093/mnras/staa903
  * Santangelo et al. (2009) Santangelo, G., Testi, L., Gregorini, L., et al. 2009, A&A, 501, 495, doi: 10.1051/0004-6361/200911688
  * Schutte & Reines (2022) Schutte, Z., & Reines, A. E. 2022, Nature, 601, 329, doi: 10.1038/s41586-021-04215-6
  * Somerville & Davé (2015) Somerville, R. S., & Davé, R. 2015, ARA&A, 53, 51, doi: 10.1146/annurev-astro-082812-140951
  * Viti et al. (2014) Viti, S., García-Burillo, S., Fuente, A., et al. 2014, A&A, 570, A28, doi: 10.1051/0004-6361/201424116
  * Volonteri et al. (2021) Volonteri, M., Habouzit, M., & Colpi, M. 2021, Nature Reviews Physics, 3, 732, doi: 10.1038/s42254-021-00364-9
  * Wilson et al. (2008) Wilson, C. D., Petitpas, G. R., Iono, D., et al. 2008, ApJS, 178, 189, doi: 10.1086/590910
  * Zhu et al. (2017) Zhu, H., Tian, W., Li, A., & Zhang, M. 2017, MNRAS, 471, 3494, doi: 10.1093/mnras/stx1580


# Zeeman dependence of the quasiparticle scattering rate and ARPES in copper
oxides and related materials

George Kastrinakis  Institute of Electronic Structure and Laser (IESL),
Foundation for Research and Technology - Hellas (FORTH), P.O. Box 1527,
Iraklio, Crete 71110, Greece∗

(29 May 2023)

###### Abstract

Within a strongly interacting Fermi liquid framework, we calculate the effects
of the Zeeman energy $\omega_{H}$ for a finite magnetic field, in a metallic
system with a van Hove peak in the density of states, located close to and
below the Fermi surface. We find that the chemical potential increases with
the square of $\omega_{H}$. We obtain a characteristic quasiparticle
scattering rate linear in the maximum of $\omega_{H}$ and temperature, both in
the normal and the d-wave superconducting state. We predict that ARPES
experiments should be able to elucidate this behavior of the scattering rate,
and in particular, the difference between spin up and down electrons.

1\. Introduction

Recently, two different angle-resolved photoemission (ARPES) experiments on
cuprates, and related materials, in finite magnetic field $H$, have appeared
ryu ; huang . ARPES experiments have the capacity to yield the quasi-particle
scattering rate, as a function of temperature $T$ and energy $\epsilon$, for
momenta anywhere close to the Fermi surface (FS). The latest experiments
include the effects of a finite $H$. High resolution ARPES experiments on the
cuprates, which demonstrated the $T,\epsilon$-linear scattering rate, started
with Valla et al. valla  in 1999 for the normal state. For the d-wave
superconducting state, the linearity of the scattering rate was pioneered in
refs. john ; kami ; yusof .

The quasi-particle scattering rate is inextricably linked to the resistivity
and the optical conductivity. In ref. gk0 ; gk1 , in 1997, for a pure Fermi
liquid (FL), with a van Hove singularity, located at energy $\epsilon_{vH}$
below the chemical potential $\mu$, we had predicted a rate linear in the
maximum of $T$ and $\epsilon$, namely for $T>x_{0}/4,\epsilon$, and for
$\epsilon>x_{0},T$, with $x_{0}=\mu-\epsilon_{vH}$. ARPES data lu  showed that
$x_{0}$ is in the range 5-30 meV (and even smaller) for several cuprates. The
$\epsilon$-linear behavior is limited above by the system bandwidth, while the
$T$-linear one is limited above by the melting of the crystal, or by some
other phase transition which may alter the value of $x_{0}$. Then this rate
yields directly a $T$-linear resistivity, and a $1/\epsilon$ optical
conductivity, for the appropriate $T,\epsilon$ ranges gk0 ; gk1 . E.g. the
scale of $x_{0}$ is fully consistent with the onset of $T$-linear resistivity.

Herein we further develop our FL model for quasiparticle scattering gk0 ; gk1
; gk2 . The extension of the model here for finite magnetic field $H$
comprises only the Zeeman energy $\omega_{H}$, and omits any dependence on
orbital effects, and on other order parameters, such as density waves of any
kind.

In section 2 we discuss the model, the chemical potential and, briefly, the
density of states (DOS) for finite $H$. These are needed in section 3, where
we discuss the quasiparticle scattering rate, and related predictions for
ARPES experiments. We conclude in section 4.

2\. Model and chemical potential for finite $\omega_{H}$

We wish to model the copper oxides. We consider the single band Hamiltonian
with a generic and non-separable potential $V_{q}$

$\displaystyle
H=\sum_{k,\sigma}(\epsilon_{k}-\mu)\;\;c^{\dagger}_{k,\sigma}c_{k,\sigma}+\frac{1}{2}\sum_{k,p,q,\sigma,\sigma^{\prime}}V_{q}\;c^{\dagger}_{k+q,\sigma}c^{\dagger}_{p-q,\sigma^{\prime}}c_{p,\sigma^{\prime}}c_{k,\sigma}\;\;.\;\;$
(1)

$c^{\dagger}_{k,\sigma}$ creates an electron of momentum $k$ and spin
$\sigma$, $\mu$ is the chemical potential, and we use the relevant band
structure

$\epsilon_{k}=-2t(\cos k_{x}+\cos k_{y})-4t^{\prime}\cos k_{x}\cos
k_{y}-2t^{\prime\prime}\\{\cos(2k_{x})+\cos(2k_{y})\\}\;\;,\;\;k_{x},k_{y}=-\pi\rightarrow\pi\;\;,\;\;$
(2)

which yields both the correct FS, and a van Hove singularity in the DOS,
located at $\epsilon=\epsilon_{vH}<\mu$, corresponding to
$k=(\pm\pi,0),(0,\pm\pi)$. The electron-electron interaction results in
pinning $\epsilon_{vH}$ slightly below $\mu$ for a broad range of dopings,
which may extend fron the underdoped to the overdoped regime in various
cuprates gk1 . Early experimental evidence from ARPES can be found in lu .

We wish to calculate $\mu_{H}=\mu(\omega_{H})$ for finite $H$. For this, we
need a self-consistent calculation of the ground state energy, taking into
account that the spin-dependent dispersion becomes

$\epsilon_{k,\sigma}=\epsilon_{k}\mp\omega_{H}\;\;,\;\;\sigma=\uparrow,\downarrow\;\;,\;\;$
(3)

and $\epsilon_{k}$ given by eq. (2). This results in the FS expanding and
shrinking for spin $\sigma=\uparrow,\downarrow$, respectively, such that the
total filling factor remains constant with $\omega_{H}>0$.

We perform a calculation which is simpler than the full Eliashberg-type one.
We work in the static limit (without any frequency dependence). We use a
variational BCS wavefunction, which allows to calculate the ground state
energy at temperature $T=0$ as

$E=\sum_{k,\sigma}(\epsilon_{k,\sigma}-\mu_{H})\;n_{k,\sigma}-\frac{1}{2}\sum_{k,q,\sigma}V_{q}\;n_{k+q,\sigma}\;n_{k,\sigma}+\frac{1}{2}\sum_{k,p,\sigma,\sigma^{\prime}}V_{q=0}\;n_{k,\sigma}\;n_{p,\sigma^{\prime}}\;\;.\;\;$
(4)

This is within the Hartree-Fock approximation. Superconducting correlations
are omitted. The occupation factors $n_{k,\sigma}$ take only the values 0 and
1 zcom . The total filling factor is $n=\sum_{k,\sigma}n_{k,\sigma}$. We use
the potential

$V_{q}=\sum_{i=1}^{4}\frac{V_{0}}{a_{0}^{2}+\xi^{2}(\vec{q}-\vec{Q}_{i})^{2}}\;\;.$
(5)

The peaks of $V_{q}$ are symetrically located at $Q_{i}=(\pm\pi,\pm\pi)$. We
typically take $\xi=2$, $a_{0}=1$, $V_{0}=6t$. $V_{q}$ may be considered as an
effective screened potential, e.g. like the potentials in the full self-
consistent solutions of the Eliashberg equations - c.f. ref. gk1  and therein.
This is a strong coupling approach, with full momentum dependence. We use a
240 $\times$ 240 discretization of the 2-d Brillouin zone. The filling $n$ is
held fixed in the calculation, which yields the self-consistent evaluation of
$n_{k,\sigma}$, $\mu_{H}$ and the respective $\sigma=\uparrow,\downarrow$
FS’s. The latter are shown in Figs. 1,2, and they are determined from the
condition

$\Xi_{k,\sigma}=0\;\;,\;\;$ (6)

where $\Xi_{k,\sigma}$ is the renormalized dispersion

$\displaystyle\Xi_{k,\sigma}=\frac{\partial E}{\partial
n_{k,\sigma}}=\epsilon_{k,\sigma}-\mu_{H}+\Sigma_{1\sigma}(k,0)\;\;,\;\;$ (7)
$\displaystyle\Sigma_{1\sigma}(k,0)=\sum_{q}\;V_{q}\;n_{k+q,\sigma}-V_{q=0}\;n\;\;.\;\;$
(8)

The condition $\Xi_{k,\sigma}>0$ yields $n_{k,\sigma}=0$. The (real) self-
energy $\Sigma_{1\sigma}(k,0)$ here is assumed to correspond to the
$\epsilon=0$ limit of the real part of the fully $\epsilon$-dependent self
energy $\Sigma_{1\sigma}(k,\epsilon)$ appearing in eq. (15) and further on.

Figure 1:  Fermi surface for $t^{\prime}=-0.2t$, $t^{\prime\prime}=0.15t$ and
$V_{q}$ of eq. (5), and filling $n=0.81$. Middle line $\omega_{H}=0$. Lower
line $\omega_{H}=0.4\;t$ spin-up, upper line $\omega_{H}=0.4\;t$ spin-down.
Also c.f. text.

Figure 2:  Fermi surface for the same parameters as above, and filling
$n=0.91$. Middle line $\omega_{H}=0$. Lower line $\omega_{H}=0.4\;t$ spin-up,
upper line $\omega_{H}=0.4\;t$ spin-down. Also c.f. text.

As can be seen in Fig. 3, we obtain a parabolic dependence of $\mu_{H}$ with
$\omega_{H}$

$\mu_{H}=\mu_{0}+a\;\omega_{H}^{2}\;\;,\;\;\omega_{H}<\omega_{*}(n)\;\;,\;\;$
(9)

where the filling dependent $\omega_{*}(n)$ is the value beyond which the FS
of the spin down electrons shrinks below the M points $(\pm\pi,0),(0,\pm\pi)$.
The constant $a$ depends on the filling factor (and other band parameters), as
shown in Table 1. Though the data shown is for a particular set of band
parameters, eq. (9) holds in general (with modified $\omega_{*}$’s etc.), as
we have seen in our numerical calculations. We also consider the case
$\mu_{H}-\mu_{0}=a\;\omega_{H}^{2}<\omega_{H}/4$. This is true for
$\omega_{H}<\omega_{0}$ with

$\omega_{0}=\frac{1}{4a}\;\;.\;\;$ (10)

In this case ($a\;\omega_{H}^{2}$) can be ignored compared to $\omega_{H}$,
and this is relevant for the subsequent discussion of the scattering rate.
Further

$\omega_{*}(n)=\omega_{0}(n)-\delta\omega(n)\;\;,\;\;\delta\omega(n)\ll\omega_{*}(n)\;\;.\;\;$
(11)

Values of $\omega_{0}$ and $\delta\omega$ also appear in Table 1.

filling n |  $\mu_{0}$ $(t)$ |  $a$ $(t^{-1})$ |  $\omega_{*}$ $(t)$ |  $\omega_{0}$ $(t)$ |  $\delta\omega=\omega_{0}-\omega_{*}$ $(t)$
---|---|---|---|---|---
$0.71$ | -1.1806 | 0.53053 | $\sim 0.4$ | 0.471 | 0.071
$0.76$ | -1.0727 | 0.40521 | $\sim 0.55$ | 0.617 | 0.067
$0.81$ | $-0.9528$ | $0.32191$ | $\sim 0.75$ | $0.777$ | 0.027
$0.835$ | $-0.8868$ | $0.28429$ | $\sim 0.85$ | $0.879$ | 0.029
$0.86$ | $-0.8195$ | $0.26061$ | $\sim 0.95$ | $0.959$ | 0.009
$0.885$ | $-0.7476$ | $0.23405$ | $\sim 1.05$ | $1.07$ | 0.02
$0.91$ | $-0.6736$ | $0.21385$ | $\sim 1.15$ | $1.17$ | 0.02

Table 1. Parameters related to $\mu_{H}$, as in eqs. (9), (10), (11), for
$t^{\prime}=-0.2t$, $t^{\prime\prime}=0.15t$ and $V_{q}$ of eq. (5), for
filling factors $n=0.71-0.91$. Also c.f. Figs. 1, 2 and text.

Figure 3:  $\mu_{H}$ $(t)$ vs. $\omega_{H}$ $(t)$ for $t^{\prime}=-0.2t$,
$t^{\prime\prime}=0.15t$ and $V_{q}$ of eq. (5). The curves bottom to top
correspond to the filling factors $n=0.71,0.76,0.81,0.835,0.86,0.885,0.91$.
The symbols come out of the self-consistent numerical calculations. The lines
are fits of eq. (9) with the parameters given in Table 1. Also c.f. text.

As shown e.g. in Fig. 4 of ref. gk1 , the DOS $N(\epsilon)$ can be
approximated as a smooth background with a single peak for $\omega_{H}=0$,
located close to and below the FS. This is known both from numerous ARPES
experiments lu  on various materials, and from several calculations (like in
gk1  and therein). For finite $\omega_{H}>0$ the DOS splits into two distinct
$\uparrow,\downarrow$ components, self adjusting in order to conserve the
total spectral weight.

The $\omega_{H}^{2}$ dependence of $\mu_{H}$ is also obtained through an
analytical calculation, wherein the DOS $N(\epsilon)$ is taken as a constant
background, with an additional high peak slightly below the FS for
$\omega_{H}=0$. Then we require conservation of the spectral weight 2
$\int_{-A}^{\mu_{0}}d\epsilon\;N(\epsilon)=\int_{-A}^{\mu_{H}}d\epsilon\;N(\epsilon+\omega_{H})+\int_{-A}^{\mu_{H}}d\epsilon\;N(\epsilon-\omega_{H})$.
Here $(-A)$ is the lower energy cut-off.

3\. Quasiparticle scattering rate

We calculate the quasiparticle scattering rate within a FL framework, as in
our previous work gk0 ; gk1 ; gk2 . Our approach is also within the framework
of conserving Baym-Kadanoff approximations bk . The derivation for the rate is
valid both in the normal and the d-wave superconducting state. In gk2  we
analysed relevant limitations for a s-wave superconducting state. Further, we
assume that any $T>0$ dependence of $\mu_{H}$ is subleading compared to the
other effects of $T>0$.

The scattering rate is a genuine FL type one, in the sense that for
$x$=max$\\{T,\epsilon\\}\rightarrow 0$, it is proportional to $x^{2}$ as in
gk0 ; gk1 ; gk2 . As a result, there is no concomitant analyticity issue (as
e.g. in Krammers-Kronig relations connecting the real and imaginary parts of
the self-energy). Only within appropriate $T$ and $\epsilon$ bounds is the
rate linear in $x$, and possibly linear in $\omega_{H}$ (as shown below).

We start from the Eliashberg equations, as e.g. written in gk2 , and follow
the same notation. The quasiparticle self-energy is
$\Sigma(k,\epsilon)=\text{Tr}\;V\;G$, where $V(q,\omega)$ is the effective
potential and $G(k,\epsilon)$ the Green’s function. We use the expression for
the scattering rate, i.e. the imaginary part of $\Sigma(k,\epsilon)$ in eq.
(18) therein

$\Sigma_{2}(k,\epsilon)\simeq\frac{1}{2}\sum_{q}V_{q}^{(1)}(0)B(k-q,\epsilon-
w_{k-q,\epsilon})\;w_{k-q,\epsilon}\;A(w_{k-q,\epsilon},\epsilon,T)\;\;.\;\;$
(12)

Here $V_{q}^{(1)}(0)=\partial\;\text{Im}V(q,\omega=0)/\partial\omega$, and we
have the thermal factor

$A(w_{k,\epsilon},\epsilon,T)=\coth\left(\frac{w_{k,\epsilon}}{2T}\right)+\tanh\left(\frac{\epsilon-
w_{k,\epsilon}}{2T}\right)\;\;,\;\;$ (13)

which satisfies

$A(w_{k,\epsilon},\epsilon,T)=0\;\;,\;\;\epsilon<w_{k,\epsilon}<0\;\;.\;\;$
(14)

The energy factor $w_{k,\epsilon}$ is given by

$w_{k,\epsilon}=\frac{\epsilon+\mu-\epsilon_{k}-\Sigma_{1}(k,\epsilon)-D_{1}(k,\epsilon)}{1-\partial_{\epsilon}\Sigma_{1}(k,\epsilon)-\partial_{\epsilon}D_{1}(k,\epsilon)}\;\;.\;\;$
(15)

$\epsilon$ is the quasiparticle energy, $\Sigma_{1}(k,\epsilon)$ is the real
part of $\Sigma(k,\epsilon)$, and $D_{1}(k,\epsilon)$ the real part of the
superconducting gap. For a small gap, we have the conditions
$D_{i}^{2}(k,\epsilon)\ll f_{i}^{2}(k,\epsilon)$ (c.f. between eqs. (10) and
(11) in gk2  for $f_{i}$ and $D_{i}$, $i=1,2$). Then we get

$B(k,\epsilon)=1-\frac{D_{1}(k,\epsilon)\;f_{1}(k,\epsilon)+D_{2}(k,\epsilon)\;f_{2}(k,\epsilon)}{f_{1}^{2}(k,\epsilon)+f_{2}^{2}(k,\epsilon)}\;\;.\;\;$
(16)

In the normal state $B(k,\epsilon)=1$.

For finite $H$ the above $k$-dependent quantities and $\mu$ become
$\omega_{H}$-dependent. Omitting the $\omega_{H}$-dependence of
$D_{1}(k,\epsilon)$ we now have

$w_{k,\epsilon,H}^{\sigma}=\frac{\epsilon+\mu_{H}-\epsilon_{k}\pm\;\omega_{H}-\Sigma_{1\sigma}(k,\epsilon)-D_{1}(k,\epsilon)}{1-\partial_{\epsilon}\Sigma_{1\sigma}(k,\epsilon)-\partial_{\epsilon}D_{1}(k,\epsilon)}\;\;.\;\;$
(17)

From our numerical calculations we see that $\Sigma_{1\sigma}(k,0)$ depends on
$\omega_{H}$, namely the part $S_{\sigma}(k)=\sum_{q}\;V_{q}\;n_{k+q,\sigma}$.
It can be very well approximated by

$S_{\sigma}(k)=S_{0}(k)+g_{\sigma}(k)\;\omega_{H}\;\;,\;\;\omega_{H}<\omega_{f}(n)\;\;,\;\;$
(18)

with the filling dependent $\omega_{f}(n)<\omega_{*}$. In table 2 we show the
FS averaged values $S_{0}=<S_{0}(k)>_{FS}$, $g_{\sigma}=<g_{\sigma}(k)>_{FS}$
of the parameters of eq. (18).

We draw the attention of the reader to the analysis of the relevant
inequalities (19)-(24) in gk2 . Repeating the analysis therein, for
$T>\epsilon,w_{k,\epsilon,H}^{\sigma}$ the scattering rate is linear in $T$

$\Sigma_{2\sigma}^{R}(k,\epsilon,H)\simeq-T\;\sum_{q}V_{q}^{(1)}(0)\;B(k-q,\epsilon-
w_{k-q,\epsilon,H}^{\sigma})\;\;.\;\;$ (19)

filling $n$ |  $S_{0}$ $(t)$ | $g_{\uparrow}$ | $g_{\downarrow}$ |  $\omega_{f}$ $(t)$
---|---|---|---|---
0.71 | 0.1867 | 0.1840 | -0.1755 | 0.4
0.76 | 0.2088 | 0.1628 | -0.1420 | 0.45
0.81 | 0.2259 | 0.1639 | -0.1379 | 0.75
0.835 | 0.2358 | 0.1615 | -0.1355 | 0.8
0.86 | 0.2463 | 0.1580 | -0.1324 | 0.8
0.885 | 0.2572 | 0.1539 | -0.1287 | 0.9
0.91 | 0.2689 | 0.1485 | -0.1247 | 1.0

Table 2. Parameters related to $S_{\sigma}$, averaged over the FS, as in eq.
(18), for $t^{\prime}=-0.2t$, $t^{\prime\prime}=0.15t$ and $V_{q}$ of eq. (5),
for filling factors $n=0.71-0.91$. Also c.f. text.

Figure 4:  $S$ $(t)$, averaged over the FS, vs. $\omega_{H}$ $(t)$ for
$t^{\prime}=-0.2t$, $t^{\prime\prime}=0.15t$ and $V_{q}$ of eq. (5). The
curves bottom to top correspond to the filling factors $n=0.71,0.81,0.91$. The
symbols come out of the self-consistent numerical calculations. The lines are
fits of eq. (18) with the parameters given in Table 2. For $n=0.71$ the linear
fit fails for smaller $\omega_{H}$ than the other 2 cases. Also c.f. text.

In the limit $\omega_{H}>T$ we obtain

$\Sigma_{2\sigma}^{R}(k,\epsilon,H)\simeq-\frac{d_{k}^{\sigma}}{2}\sum_{q}V_{q}^{(1)}(0)\;B(k-q,\epsilon-
w_{k-q,\epsilon,H}^{\sigma})\;w_{k-q,\epsilon,H}^{\sigma}\;\;.\;\;$ (20)

Here $d_{k}^{\sigma}$=1+$<$tanh[$(\epsilon-
w_{k-q,\epsilon,H}^{\sigma})/(2T)$]$>$ is a factor between 1 and 2. Also note
that

$\partial_{\epsilon}\Sigma_{1\sigma}(k,\epsilon)=-b_{k}<0\;\;,\;\;\epsilon<\epsilon_{o}=O(t)\;\;.\;\;$
(21)

Then we obtain

$\Sigma_{2\sigma}^{R}(k,\epsilon,H)\simeq-\frac{d_{k}^{\sigma}}{2}\sum_{q}V_{q}^{(1)}(0)\;B(k-q,\epsilon-
w_{k-q,\epsilon,H}^{\sigma})\;\frac{\epsilon+\mu_{H}-\epsilon_{k-q}\pm\;\omega_{H}-\Sigma_{1\sigma}(k-q,\epsilon)}{1+b_{k-q}}\;\;.\;\;$
(22)

There is interference between $\mu_{H}$, the pure Zeeman term
$\pm\;\omega_{H}$, and the term $g_{\sigma}\;\omega_{H}$. It should be
possible to probe the rates $\Sigma_{2\sigma}^{R}(k,\epsilon,H)$ by performing
ARPES experiments in adequately strong field $H$. In particular, it should
also be possible to measure the pure Zeeman dependence, with the assumption
that $d_{k}^{\sigma}=d_{k}$ has a negligible dependence on $\omega_{H}$. The
appropriate quantity is the difference

$\Sigma_{2\uparrow}^{R}(k,\epsilon,H)-\Sigma_{2\downarrow}^{R}(k,\epsilon,H)=-d_{k}\;\omega_{H}\;\sum_{q}\frac{V_{q}^{(1)}(0)}{1+b_{k-q}}\;B(k-q,\epsilon-
w_{k-q,\epsilon})\;l_{k-q}\;\;,\;\;$ (23)

with

$l_{k}=1+\frac{1}{2}\\{g_{\uparrow}(k)-g_{\downarrow}(k)\\}\;\;.\;\;$ (24)

In particular, in the normal state this formula becomes

$\Sigma_{2\uparrow}^{R}(k,\epsilon,H)-\Sigma_{2\downarrow}^{R}(k,\epsilon,H)=-y_{k}\;\omega_{H}\;\;,\;\;\;$
(25)

with

$y_{k}=d_{k}\;\sum_{q}\frac{V_{q}^{(1)}(0)}{1+b_{k-q}}\;l_{k-q}\;\;.\;\;$ (26)

This is to be compared with the total prefactor for the case
$T>\epsilon,\omega_{H}$ in $\Sigma_{2}^{R}(k,\epsilon)=-x_{k}\;T$, c.f. eq.
(19),

$x_{k}=\sum_{q}V_{q}^{(1)}(0)\;\;.\;\;$ (27)

Here the two prefactors $y_{k}$ and $x_{k}$ are similar in magnitude. This
should be visible in ARPES experiments.

Further, we examine the quantity

$s_{\sigma}(k,H)=\sum_{q}V_{q}^{(1)}(0)\;\Xi_{k-q,\sigma}\;\;,\;\;$ (28)

which satisfies $\Sigma_{2\sigma}^{R}(k,\epsilon=0,H)\propto s_{\sigma}(k,H)$.
Assuming that the $\omega$-dependent $V(q,\omega)$ (compatible with $V_{q}$ in
eq. (5) ) is

$V(q,\omega)=\sum_{i=1}^{4}\frac{V_{0}}{a_{0}^{2}+\xi^{2}(\vec{q}-\vec{Q}_{i})^{2}-i\;\omega/\omega_{F}}\;\;,\;\;$
(29)

we have

$V_{q}^{(1)}(0)=\frac{V_{0}}{\omega_{F}}\;\sum_{i=1}^{4}\frac{1}{\big{[}a_{0}^{2}+\xi^{2}(\vec{q}-\vec{Q}_{i})^{2}\big{]}^{2}}\;\;.\;\;$
(30)

The small fermionic energy is typically $\omega_{F}=10-40$ meV gk1 . Here we
take $\omega_{F}=t/10$. For e.g. $t=250$ meV, this yields $\omega_{F}=25$ meV.

The numerical calculations yield

$s_{\sigma}(k,H)=s_{0}(k)+f_{\sigma}(k)\;\omega_{H}\;\;,\;\;\omega_{H}<\omega_{1}\;\;,\;\;$
(31)

typically with $\omega_{1}\leq\omega_{*}(n)$, in accordance with eq. (18). The
FS averaged coefficients $s_{\sigma}$ and $f_{\sigma}$ are shown in Table 3.
Note that $f_{\uparrow}<0$ and $f_{\downarrow}>0$.

filling $n$ | $s_{0}$ | $f_{\uparrow}$ | $f_{\downarrow}$
---|---|---|---
0.71 | 1.8021 | -1.8039 | 1.1347
0.76 | 1.6147 | -1.7839 | 1.2399
0.81 | 1.4095 | -1.7604 | 1.3055
0.86 | 1.1713 | -1.7107 | 1.3734
0.91 | 0.9224 | -1.7086 | 1.4195

Table 3. Parameters related to $s_{\sigma}$, averaged over the FS, as in eq.
(28), for $t^{\prime}=-0.2t$, $t^{\prime\prime}=0.15t$ and $V_{q}^{(1)}(0)$ of
eq. (30), with $\omega_{F}=t/10$ and for filling factors $n=0.71-0.91$.
$f_{\uparrow}<0$ and $f_{\downarrow}>0$. Also c.f. text.

Figure 5:  $s_{\sigma}$ $(t)$, averaged over the FS, vs. $\omega_{H}$ $(t)$ in
eq. (28) for $t^{\prime}=-0.2t$, $t^{\prime\prime}=0.15t$ and $V_{q}^{(1)}(0)$
of eq. (30). The curves top to bottom correspond to the filling factors
$n=0.71,0.81,0.91$. The symbols come out of the self-consistent numerical
calculations. The lines are fits of eq. (31) with the parameters given in
Table 3. For $n=0.71$ the linear fit fails for smaller $\omega_{H}$ than the
other 2 cases. Also c.f. text.

4\. Overview

Based on our earlier work, we calculate the scattering rate in a finite
magnetic field, taking only the Zeeman energy $\omega_{H}$ dependence into
account, within a purely FL framework. The DOS has a strong van Hove peak.

We make specific predictions for the quasi-particle scattering rate, which can
be probed by ARPES experiments. We find a characteristic linear in max
$\\{T,\omega_{H}\\}$ dependence. We note that there is an asymmetry between up
and down spin electrons for finite $\omega_{H}$, as is evident from figures 4
and 5, and from tables 2 and 3. This asymmetry is also filling factor $n$
dependent, increasing for lower values of $n$.

The related compounds BaFe2As2 bafe  and CeCoIn5 ceco  have shown a
combination of $T$-linear resistivity, a van Hove singularity located close to
the FS, and a $T,\epsilon$-linear scattering rate (modified possibly by
structural phase transitions). Our results apply in these compounds as well.

∗ e-mail<EMAIL_ADDRESS><EMAIL_ADDRESS>
## References

  * (1) S.H. Ryu, G. Reichenbach, C.M. Jozwiak, A. Bostwick, P. Richter, T. Seyller, and E. Rotenberg, arXiv:2304.07431.
  * (2) J. Huang, Z. Yue, A. Baydin, H. Zhu, H. Nojiri, J. Kono, Y. He, and M. Yi, arXiv:2305.02625.
  * (3) T. Valla, A. V. Fedorov, P. D. Johnson, B. O. Wells, S. L. Hulbert, Q. Li, G. D. Gu, and N. Koshizuka, Science 285, 2110 (1999).
  * (4) P. D. Johnson, T. Valla, A. V. Fedorov, Z. Yusof, B. O. Wells, Q. Li, A. R. Moodenbaugh, G. D. Gu, N. Koshizuka, C. Kendziora, S. Jian, and D. G. Hink, Phys. Rev. Lett. 87, 177007 (2001).
  * (5) A. Kaminski, M. Randeria, J. C. Campuzano, M. R. Norman, H. Fretwell, J. Mesot, T. Sato, T. Takahashi, and K. Kadowaki, Phys. Rev. Lett. 86, 1070 (2001).
  * (6) Z. M. Yusof, B. O. Wells, T. Valla, A. V. Fedorov, P. D. Johnson, Q. Li, C. Kendziora, S. Jian, and D. G. Hinks, Phys. Rev. Lett. 88, 167006 (2002).
  * (7) G. Kastrinakis, cond-mat/9703141.
  * (8) G. Kastrinakis, Physica C 340, 119 (2000).
  * (9) D. H. Lu, M. Schmidt, T. R. Cummins, S. Schuppler, F. Lichtenberg, and J. G. Bednorz, Phys. Rev. Lett. 76, 4845 (1996).
  * (10) G. Kastrinakis, Phys. Rev. B 71, 014520 (2005).
  * (11) As a result of $n_{k,\sigma}=0$ or 1 only, it is impossible to have a quasiparticle weight $z_{k}<1$ around the FS. However, since we are not concerned about such fine detail here, this is acceptable.
  * (12) G. Baym, Phys. Rev. 127, 1391 (1962); G. Baym and L. P. Kadanoff, 124, 287 (1961).
  * (13) M. Yi, D. H. Lu, J. G. Analytis, J.-H. Chu, S.-K. Mo, R.-H. He, R. G. Moore, X. J. Zhou, G. F. Chen, J. L. Luo, N. L. Wang, Z. Hussain, D. J. Singh, I. R. Fisher, and Z.-X. Shen, Phys. Rev. B 80, 024515 (2009).
  * (14) Q. Y. Chen, D. F. Xu, X. H. Niu, J. Jiang, R. Peng, H. C. Xu, C. H. P. Wen, Z. F. Ding, K. Huang, L. Shu, Y. J. Zhang, H. Lee, V. N. Strocov, M. Shi, F. Bisti, T. Schmitt, Y. B. Huang, P. Dudin, X. C. Lai, S. Kirchner, H. Q. Yuan, and D. L. Feng, Phys. Rev. B 96, 045107 (2017).