# Discovery of multiple p-mode pulsation frequencies in the roAp star, HD
86181

Fangfei Shi1,2, Donald W. Kurtz3,4, Daniel L. Holdsworth4, Hideyuki Saio5,
Margarida S. Cunha6, Huawei Zhang1,2, Jianning Fu7, G. Handler8
1Department of Astronomy, School of Physics, Peking University, Beijing
100871, P. R. China
2Kavli Institute for Astronomy and Astrophysics, Peking University, Beijing
100871, P. R. China
3Centre for Space Research, Physics Department, North-West University,
Mahikeng 2745, South Africa
4Jeremiah Horrocks Institute, University of Central Lancashire, Preston PR1
2HE, UK
5Astronomical Institute, Graduate School of Science, Tohoku University, Sendai
980-8578, Japan
6Instituto de Astrofísica e Ciências do Espaço, Universidade do Porto, CAUP,
Rua das Estrelas, PT4150-762 Porto, Portugal
7Department of Astronomy, Beijing Normal University, Beijing 100875, P. R.
China
8Nicolaus Copernicus Astronomical Center, Polish Academy of Sciences, ul.
Bartycka 18, 00-716, Warsaw, Poland

(Accepted XXX. Received YYY; in original form ZZZ)

###### Abstract

We report the frequency analysis of a known roAp star, HD 86181 (TIC
469246567), with new inferences from TESS data. We derive the rotation
frequency to be $\nu_{\rm rot}=0.48753\pm 0.00001$ d-1. The pulsation
frequency spectrum is rich, consisting of two doublets and one quintuplet,
which we interpret to be oblique pulsation multiplets from consecutive, high-
overtone dipole, quadrupole and dipole modes. The central frequency of the
quintuplet is $232.7701$ d-1 (2.694 mHz). The phases of the sidelobes, the
pulsation phase modulation, and a spherical harmonic decomposition all show
that the quadrupole mode is distorted. Following the oblique pulsator model,
we calculate the rotation inclination, $i$, and magnetic obliquity, $\beta$,
of this star, which provide detailed information about the pulsation geometry.
The $i$ and $\beta$ derived from the best fit of the pulsation amplitude and
phase modulation to a theoretical model, including the magnetic field effect,
slightly differ from those calculated for a pure quadrupole, indicating the
contributions from $\ell=4,6,8,...$ are small. Non-adiabatic models with
different envelope convection conditions and physics configurations were
considered for this star. It is shown that models with envelope convection
almost fully suppressed can explain the excitation at the observed pulsation
frequencies.

###### keywords:

stars: oscillations – stars: variables – stars: individual HD 86181 (TIC
469246567; V437 Car) – star: chemically peculiar – techniques: photometric –
asteroseismology

††pubyear: 2019††pagerange: Discovery of multiple p-mode pulsation frequencies
in the roAp star, HD 86181–LABEL:lastpage

##  1 Introduction

The Ap (chemically peculiar A-type) stars have non-uniform distributions of
chemical abundances on their surfaces and strong magnetic fields. These
magnetic fields suppress surface convection that then leads to element
stratification. For some heavy elements, such as Eu, Sr and Si, the radiation
pressure can lift them up to the surface against gravity leading to many
absorption features. These elemental overabundances occur in spots, making Ap
stars obliquely rotating variable stars of a class known as $\alpha^{2}$ CVn
stars (Pyper, 1969).

Some cool Ap stars exhibit high-overtone, low-degree pressure pulsation modes
with periods between 4.7 and 24 min (frequencies in the range $55.8-300$ d-1;
$0.6-3.5$ mHz (Holdsworth et al., 2021)) and photometric amplitudes up to
0.018 mag in Johnson $B$ (Cunha et al., 2019; Kochukhov, 2009; Smalley et al.,
2015). They are called rapidly oscillating Ap (roAp) stars. Some of these
stars show both rotation features with periods of days to decades, and
pulsation features in their light curves.

Stibbs (1950) developed the oblique rotator model of the Ap stars, which
accounts for the magnetic, spectral, and light variations observed in Ap
stars. Following this model, Kurtz (1982) introduced the oblique pulsator
model, which was generalized with the effects of both the magnetic field and
rotation taken into account (Kurtz, 1982; Dziembowski & Goode, 1985;
Shibahashi & Takata, 1993; Takata & Shibahashi, 1994, 1995; Saio & Gautschy,
2004; Bigot & Dziembowski, 2002; Bigot & Kurtz, 2011). According to this
model, the pulsation axis is misaligned with the rotation axis, and generally
closely aligned to the magnetic axis. When the star rotates, the viewing
aspect of the pulsation modes varies along the line of sight, leading to
apparent amplitude and phase modulation. This modulation can provide
information on the geometry of observed pulsations, hence mode identification,
which is necessary for asteroseismic inference with forward modelling.

Since the first roAp stars were discovered by Kurtz (1982), 88 roAp stars have
been found (Smalley et al., 2015; Hey et al., 2019; Cunha et al., 2019;
Balona, Holdsworth & Cunha, 2019; Holdsworth et al., 2021). Asteroseismology
is a useful method to diagnose stellar structure and interior physics from the
evidence of surface pulsations (Cunha, Fernandes & Monteiro, 2003). Progress
of this research for roAp stars has been hindered by the relatively small
number of known stars, and because their rapid pulsation requires dedicated
observations and high accuracy to detect the small pulsation amplitudes (Hey
et al., 2019; Cunha et al., 2019; Balona, Holdsworth & Cunha, 2019).

The space telescopes Kepler and the TESS (Transiting Exoplanet Survey
Satellite) provide an opportunity to detect oscillations well below the
amplitude threshold of ground-based observations. Both Kepler and TESS have
short cadence (2 min for TESS and 58.89 s for Kepler) observations, but Kepler
only observed 512 stars in this mode during each observing ‘quarter’. However,
the standard long cadence sampling frequency of the Kepler 30-min observations
is generally too low for studying the pulsation of roAp stars in detail.

Murphy, Shibahashi & Kurtz (2013) showed that the Nyquist ambiguity in the LC
data can be resolved as a result of the Barycentric corrections applied to
Kepler time stamps, and Hey et al. (2019) discovered 6 roAp candidates through
this method. Compared to the Kepler 58.89-s observations, TESS is observing
many more stars with 2-min observations with sufficiently long time bases to
detect pulsations. Up to now, 21 new roAp stars have been found from just TESS
sectors 1 to 13 (Cunha et al., 2019; Balona, Holdsworth & Cunha, 2019;
Holdsworth et al., 2021).

Before the TESS observations of our target, HD 86181, Kurtz & Martinez (1994)
discovered it to be a roAp star from 4.85 hr of ground-based data. They
reported the star to have a pulsation period of 6.2 min and with an amplitude
of 0.35 mmag through a Johnson $B$ filter. That period corresponds to a
frequency of 2.688 mHz, or 232.26 d-1. No further detailed studies of the
pulsations in HD 86181 have been published.

Parameters for this star are listed in Table 1. The effective temperature was
estimated using the Str$\ddot{o}$mgren photometric indices extracted from the
catalogue of Hauck & Mermilliod (1998) and the calibrations in the TEMPLOGG
code (Rogers 1995) which were developed based on the work of Moon & Dworetsky
(1985) and Napiwotzki, Schoenberner & Wenske (1993). Since no convincing
uncertainty is given by this method, we indicate, instead, a range of values
of $T_{\rm eff}$ published in the literatures.

The luminosity was calculated through the relation $-2.5\log
L=M_{G}+BC_{G}(\rm T_{eff})-M_{bol,\odot}$, where $BC_{G}(\rm T_{eff})$ is a
temperature dependent bolometric correction defined in Andrae et al. (2018),
and the uncertainty of BC (Bolometric Correction) is 0.13, based on a
comparison with Ap data that is described in some detail in Cunha et al.
(2019). While the uncertainty derived in Cunha et al. (2019) was based on a
comparison of Ap-star measurements with the empirical BCV calibration by
Flower (1996) and, thus, the consistency of using it with the BGG values
derived from the calibration of Andrae et al. (2018) may be questionable, it
provides a more conservative result than the uncertainty derived from Andrae
et al. (2018), which does not account for the stars’ peculiarities. The
extinction in the $G$ band used to calculate $M_{G}$ here was from Anders et
al. (2019), and the uncertainty is 0.2, which is the value indicated in the
Figure 20 in (Anders et al., 2019). The parallax was from Gaia eDR3 (Gaia
Collaboration, 2020). $M_{\rm bol,\odot}$ adopted is 4.74 as defined by IAU
Resolution 2015
B2111https://www.iau.org/static/resolutions/IAU2015_English.pdf.

Table 1: Parameters of HD 86181.

Apparent $G$ magnitude | $9.341\pm 0.003$ | Gaia Collaboration (2020)
---|---|---
Extinction in $G$ band | $0.1\pm 0.2$ | Anders et al. (2019)
Spectral type | F0 Sr | Renson & Manfroid (2009)
Parallax (mas) | $4.15\pm 0.01$ | Gaia Collaboration et al. (2020)
Distance (pc) | $241.0\pm 0.6$ | derived from parallax
$b-y$ | 0.175 | Perry (1991)
$m_{1}$ | 0.245 | Perry (1991)
$c_{1}$ | 0.702 | Perry (1991)
Hβ | 2.804 | Perry (1991)
$T_{\rm eff}$(K) |  $7750$; [7240-7910] |  This work∗; Literature+
Luminosity (L⊙) | $8.8\pm 1.9$ | Andrae et al. (2018)
Mean longitudinal | $536\pm 75$ | Bagnulo et al. (2015)
magnetic field (G) |  |

  * •

∗ based on Rogers 1995

  * •

+ Trifonov et al. (2020), Anders et al. (2019), Mathys, Kharchenko & Hubrig (1996)

##  2 TESS observations

HD 86181 was observed by TESS in sectors 9 and 10 in 2-min cadence. The data
have a time span of 51.76 d with a centre point in time of $t_{0}={\rm
BJD}~{}2458569.80077$, and comprise 33832 data points after some outliers were
removed. The standard PDC SAP (pre-search data conditioning simple aperture
photometry) fluxes provided by MAST (Mikulski Archive for Space Telescopes)
were used and normalised by dividing by the median flux separately for each
sector. Relative magnitudes were then calculated from the processed fluxes,
giving the light curve shown in the top panel of Fig. 1.

There are obvious rotational variations from spots, as is typical of the
magnetic Ap stars. Within the oblique rotator model, the double wave nature of
the rotational variations suggests that two principal spots with enhanced
brightness on the stellar surface are seen. The high frequency pulsation
cannot be seen in this figure at this resolution.

##  3 Frequency analysis

###  3.1 Rotation frequency analysis

Before we conducted a detailed analysis of the rotation frequency of HD 86181,
we first measured the rotation frequency with a coarse Discrete Fourier
Transform (DFT; Kurtz, 1985) such that we could bin the data every one
rotation cycle. This allowed us to assess the instrumental variation, which we
subsequently fit with a polynomial and removed from the original light curve.
We then calculated a DFT with a finer frequency grid, as shown in Fig. 2, to
measure the stellar rotation frequency. The low frequencies dominate in the
spectrum, so we zoom in to both the low frequency range (second panel) and
high frequency range (third panel). From the amplitude spectrum at low
frequency, the rotational harmonics are clearly seen. Although the highest
peak is at a frequency of 0.97 d-1, considering the phase plot, we derive the
rotation frequency to be around 0.48 d-1. Because the variation is a double
wave, the second harmonic has the highest amplitude.

A linear least-squares fit was calculated to find the best amplitudes and
phases of the rotation frequency and its 4 visible harmonics, and then a non-
linear least-squares fit to get optimized results. The rotational frequency is
derived to be $\nu_{\rm rot}=0.48753\pm 0.00001$ d-1 ($P_{\rm rot}=2.05116\pm
0.00004$ d) by dividing the frequency of the highest amplitude second harmonic
by two, which has better signal-to-noise ratio. Besides the rotation frequency
and its harmonics, there are still some signals left in the low frequency
range, probably because instrumental variation has not been removed
completely. These signals were removed prior to the non-linear least square
fits for better estimates of the uncertainties. The uncertainties were derived
following Montgomery & O′donoghue (1999). The rotation period is short among
the known roAp stars, after HD 43226 (Cunha et al. 2019), HD 216641 (Cunha et
al. 2019), and HD 6532 (Kurtz et al. 1996a, Kurtz et al. 1996b), which have
similar rotation periods of $P_{\rm rot}=1.71441$ d, $P_{\rm rot}=1.876660$ d,
and $P_{\rm rot}=1.944973$ d, respectively.

Figure 1: Top: The light curve of HD 86181 showing the rotational variations.
Bottom: Phase folded light curve of HD 86181, folded on the rotation period of
2.05116 d; two rotation cycles are shown for clarity. The data are from TESS
sectors 9 and 10. The time zero-point, BJD 2458569.26128, is the time of
pulsation maximum. The phases are binned every 0.001 phase bin.

Figure 2: The frequency spectrum of HD 86181. Top: The amplitude spectrum of
the S9–10 data out to 300 d-1. The rotational frequencies at low frequency
dominate. The pulsation frequencies centred on 232.2 d-1 are difficult to see
at this scale. Second: The low frequency rotational harmonics. Third: the
pulsation frequencies for the high-pass filtered data. Bottom: The frequency
spectrum after the frequencies in Table 2 have been removed. The red
horizontal lines are 4 times of noise level. The top x-axis is the
corresponding frequency in $\mu$Hz.

###  3.2 The pulsations

To study the pulsations, a high-pass filter was used to remove the rotational
light variations, any remaining instrumental artefacts and other low
frequencies. The high-pass filter was a simple consecutive pre-whitening of
low frequency peaks extracted by Fourier analysis until the noise level was
reached in the frequency range $0-6$ d-1. The third panel in Fig. 2 shows the
amplitude spectrum for the high-pass filtered data around the high-frequency
variability. By inspection it can be seen that there is a central quintuplet
and two doublets, one at higher and another at lower frequency than the
quintuplet. After removing these three groups of frequencies, five singlets
still remain (see the bottom panel of Fig. 2). However, their frequencies are
similar to the quintuplet and two doublets within the uncertainties. These may
be caused by amplitude or frequency modulation over the time span of the data
set, 51.76 d.

To test this, we removed the doublets and singlets from the light curve and
fitted $\nu_{1}$ to sections of the data that are exactly one rotation cycle
long and calculated the amplitude and phase. Fig. 3 shows there is amplitude
and phase variability with time. By choosing exactly one rotation length of
data, the amplitude and phase variations due to oblique pulsation were
smoothed.

If the frequency were stable, there would be no phase variations. As the data
were fitted with the function $\Delta m=A\cos(\nu(t-t_{0})+\phi)$, the
frequency and phase terms are inextricably intertwined (see the section 5.3.2
in Holdsworth et al., 2014), thus a change in one can be interpreted as a
change in the other. Therefore, although we show a change in the phase in Fig.
3, the change could be in the frequency. Such variability is common in roAp
stars studied with high precision data (Holdsworth, 2021).

Figure 3: The pulsation amplitude and phase variations of HD 86181 for the
dominant quadrupole mode. Top: pulsation amplitude variations as a function of
time. Bottom: pulsation phase variations as a function of time.

As in the analysis of rotation frequency, linear and non-linear least squares
fits were used to get optimised results of frequencies, amplitudes and phases.
The non-linear least squares fit results are shown in Table 2. Within the
uncertainties, the sidelobes of the quintuplet are exactly split by the
rotation frequency. In addition to the quintuplet, there are two doublets that
are split by 2$\nu_{\rm rot}$; these are the sidelobes of two dipole pulsation
frequencies that are labeled as $\nu_{2}$ and $\nu_{3}$. For a pure dipole or
quadrupole pulsation, the oblique pulsator model requires that the sidelobes
are split by exactly the rotation frequency of the star, and that the phases
of all components are equal at the time of pulsation maximum. To test this,
the frequency of the quintuplet sidelobes were fixed to be equally spaced by
the rotation frequency, and the zero-point in time was chosen such that the
phases of the first pair of sidelobes are the same, then a linear least
squares fit was applied to the data with the results show in Table 3. The
phases of the quintuplet sidelobes are not equal within the uncertainties,
which indicates this star pulsates in a distorted quadrupole mode.

Table 2: A non-linear least squares fit of the frequency multiplets for HD 86181. The zero point for the phases is $t_{0}={\rm BJD}~{}2458569.26128$. | frequency | amplitude | phase
---|---|---|---
| d-1 | mmag | radians
|  | $\pm 0.007$ |
$\nu_{rot}$ | $0.48765\pm 0.00003$ | 2.970 | $5.829\pm 0.003$
$2\nu_{rot}$ | $0.97506\pm 0.00001$ | 7.296 | $2.776\pm 0.001$
$3\nu_{rot}$ | $1.46233\pm 0.00013$ | 0.732 | $0.312\pm 0.013$
$4\nu_{rot}$ | $1.95043\pm 0.00013$ | 0.761 | $6.177\pm 0.013$
$6\nu_{rot}$ | $2.92585\pm 0.00050$ | 0.190 | $0.215\pm 0.049$
$\nu_{2}-\nu_{\rm rot}$ | $229.6162\pm 0.0012$ | 0.059 | $0.28\pm 0.17$
$\nu_{2}+\nu_{\rm rot}$ | $230.5897\pm 0.0014$ | 0.050 | $0.49\pm 0.20$
$\nu_{1}-2\nu_{\rm rot}$ | $231.7947\pm 0.0008$ | 0.091 | $6.00\pm 0.11$
$\nu_{1}-\nu_{\rm rot}$ | $232.2853\pm 0.0013$ | 0.055 | $6.27\pm 0.18$
$\nu_{1}$ | $232.7701\pm 0.0003$ | 0.273 | $6.24\pm 0.04$
$\nu_{1}+\nu_{\rm rot}$ | $233.2587\pm 0.0011$ | 0.062 | $6.17\pm 0.16$
$\nu_{1}+2\nu_{\rm rot}$ | $233.7438\pm 0.0008$ | 0.080 | $0.14\pm 0.12$
$\nu_{3}-\nu_{\rm rot}$ | $235.2495\pm 0.0010$ | 0.071 | $6.11\pm 0.14$
$\nu_{3}+\nu_{\rm rot}$ | $236.2261\pm 0.0012$ | 0.063 | $5.94\pm 0.16$
Table 3: A least squares fit of the frequency multiplets for HD 86181, where the frequency splitting of the rotational sidelobes has been forced to be exactly the rotation frequency. The zero point for the phases, $t_{0}={\rm BJD}~{}2458569.26128$, has been chosen to be a time when the first two orbital sidelobes of the quintuplet have equal phase. | frequency | amplitude | phase
---|---|---|---
| d-1 | mmag | radians
|  | $\pm 0.007$ |
$\nu_{2}-\nu_{\rm rot}$ | 229.6162 | 0.058 | $0.25\pm 0.11$
$\nu_{2}+\nu_{\rm rot}$ | 230.5913 | 0.049 | $0.40\pm 0.14$
$\nu_{1}-2\nu_{\rm rot}$ | 231.7950 | 0.091 | $-0.26\pm 0.07$
$\nu_{1}-\nu_{\rm rot}$ | 232.2826 | 0.053 | $-0.13\pm 0.12$
$\nu_{1}$ | 232.7701 | 0.273 | $-0.05\pm 0.02$
$\nu_{1}+\nu_{\rm rot}$ | 233.2576 | 0.061 | $-0.13\pm 0.11$
$\nu_{1}+2\nu_{\rm rot}$ | 233.7452 | 0.080 | $0.14\pm 0.08$
$\nu_{3}-\nu_{\rm rot}$ | 235.2494 | 0.071 | $-0.13\pm 0.09$
$\nu_{3}+\nu_{\rm rot}$ | 236.2245 | 0.063 | $-0.11\pm 0.11$

We also investigated the impact of the spots on the pulsations. From the
second panel of Fig. 1, the rotational variation caused by the spots amounts
to 20 ppt peak-to-peak. We therefore expect the modulation of the pulsation
caused by spots to be also a factor of 0.02 of the pulsation amplitude, which
is down to $\mu$mag, much below the noise level. So the effect of spots on the
pulsation amplitude is negligible.

Finally, harmonics of the pulsation frequencies were also searched for beyond
the Nyquist frequency, $\nu_{Ny}=359.804$ d-1. Only three similar alias groups
centred at $2\nu_{Ny}-\nu_{1}$, $2\nu_{Ny}-\nu_{2}$ and $2\nu_{Ny}-\nu_{3}$
were found, with no evidence of harmonics of pulsation frequencies.

###  3.3 Pulsation amplitude and phase modulation

To study the rotation modulation of the pulsation amplitudes and phases, the
light curve was divided into 217 segments each containing 50 pulsation cycles,
thus each segment had a time span of 0.21d, or 0.1 of a rotation cycle. Linear
least-squares fitting was applied to these segments at fixed frequency,
$\nu_{1}=232.7701$ d-1, to calculate the pulsation amplitude and phase as a
function of rotation phase. Fig. 4 shows these modulations along with the
rotation light variations for comparison.

Figure 4: The pulsation modulation for pulsation frequency $\nu_{1}$ of HD
86181. Top: The phase folded rotation light curve. Middle: pulsation amplitude
variations as a function of rotation phase. Amplitude points with $1\sigma$
errors greater than 0.12 mmag are not plotted here. Bottom: pulsation phase
variations as a function of rotation phase. Phase points with $1\sigma$ errors
greater than 0.8 rad are not plotted here. The red lines are theoretical
amplitude modulation modelled following Kurtz (1992) with the components from
Table 4. The blue line was calculated based on an oblique pure quadruple mode
(see section 4). Two rotation cycles are shown. The time zero-point is
$t_{0}={\rm BJD}~{}2458569.26128$.

The maxima of the pulsation amplitude depend on the aspect of the pulsation
axis, while the light extrema depend on the spots. The difference between the
occurrences of the extrema of the pulsation amplitude and the rotational light
variations indicates the position of spots relative to the pulsation axis. In
many Ap stars the surface positions where spots form – particularly for the
rare earth elements – is related to the magnetic field. In the case that the
spots are centred on the pulsation axis which is also fixed close to the
magnetic axis, the rotation phase of pulsation maximum coincides with, or is
near to, the rotation phase of the light extrema. As Handler et al. (2006)
showed for HD 99563, the maximum of pulsation amplitude coincides with the
maximum of rotation light in red filters, and the minimum in blue filters. The
antiphase variations in blue and red filters are related to the flux
redistribution from UV to optical caused by line blocking (Leckrone, 1973).

For HD 86181, pulsation amplitude maximum coincides with the secondary maximum
of the light curve, and after half a cycle, the secondary maximum of pulsation
amplitude coincides with the maximum of the light curve. For a pure quadrupole
pulsator, the intrinsic pulsation amplitude peaks at both pulsation poles and
at the equator. The pulsation maximum at the poles is twice that at the
equator, but with inverse phase. We assume the maximum of pulsation amplitude
is generated at the pole, while the secondary maximum by equator. This
assumption is verified with the the oblique pulsator model below.

At rotation phase 0, which we chose to be the time of pulsation maximum for
the quadrupole mode, we see that the spots show the secondary rotational light
maximum. In contrast, for another roAp star with a quadrupole mode, KIC
10685175, the maximum of the pulsation amplitude coincides with the minimum of
the rotational light (Shi et al., 2020) (Fig. 5).

Figure 5: Same as Fig. 4 for KIC 10685175. The time zero-point is $t_{0}={\rm
BJD}~{}2458711.21931$.

The pulsation phase as a function of rotation does not show a $\pi$-rad phase
reversal expected at the times of amplitude minima as would be the case for an
undistorted mode, although the pulsation phase shows bumps at those times.
This then argues for a distorted quadrupole mode, and also is similar to what
is observed in other roAp stars with well-studied quadrupole modes
(Holdsworth, Saio & Kurtz, 2019; Holdsworth et al., 2018b, c, a, 2014; Kurtz
et al., 1996b; Holdsworth et al., 2016).

We also checked the pulsation amplitude and phase modulations of the two
central frequencies ($\nu_{2}=230.1038$ d-1 and $\nu_{3}=235.7370\,$d-1) of
the dipole modes, as seen in Fig. 6. However, because of the low amplitudes,
the modulation curves are quite scattered, especially the phase modulation
curve. The two dipole modes show similar behaviour: the pulsation amplitude
reaches primary and secondary maximum at rotation phases 0 and 0.5,
respectively, the same as the quadrupole mode. The pulsation phase variations
have large errors, hence $\pi$-rad pulsation phase changes at rotation phases
0.25 and 0.75 – typical behaviour for dipole modes – are neither ruled out,
nor supported by the plots in Fig. 6.

Figure 6: Top panel: The pulsation amplitude (left) and phase (right)
modulation of the dipole central frequency $\nu_{2}=230.1038$ d-1. Bottom
panel: The pulsation amplitude (left) and phase (right) modulation of the
dipole central frequency $\nu_{3}=235.7370$ d-1. Phase points with $1\sigma$
errors greater than 1.0 rad are not plotted here. The red lines are
theoretical amplitude modulation modelled following Kurtz (1992) with the
components from Table 4. The time zero-point is $t_{0}={\rm
BJD}~{}2458569.26128$.

##  4 Oblique pulsator model

The oblique pulsator model describes the pulsation pattern of an oblique
pulsator and only considers the surface geometry of non-radial pulsation
modes. However, some spectroscopic observations (e.g. Kochukhov 2006;
Freyhammer et al. 2009) and simulations (e.g. Khomenko & Kochukhov 2009) have
shown that properties of pulsations change rapidly with height in the stellar
atmosphere and modes are substantially distorted by the magnetic field. Sousa
& Cunha (2011) and Quitral-Manosalva, Cunha & Kochukhov (2018a) have also
studied this extensively theoretically.

Recently, TESS observations of HD 6532 and HD 80316 (Holdsworth et al., 2021)
have shown that there are changes in multiplet structure comparing to the
former ground-based $B$ observations. The TESS filter is broad-band white-to-
red, which probes to a different depth in the stellar atmosphere than the $B$
filter. These new observations show the complexity of roAp pulsations and
importance of the vertical dimension. Nevertheless, the oblique pulsator model
still allows us a simple first look at the geometry of the pulsation modes.

For a normal quadrupole pulsator, the ratio of the sidelobes to the central
peak can be calculated with eqns 8 and 10 from Kurtz (1992):

$\frac{A_{+1}+A_{-1}}{A_{0}}=\frac{12\sin\beta\cos\beta\sin i\cos
i}{(3\cos^{2}\beta-1)(3\sin^{2}i-1))}$ (1)

and

$\frac{A_{+2}+A_{-2}}{A_{0}}=\frac{3\sin^{2}\beta\sin^{2}i}{(3\cos^{2}\beta-1)(3\sin^{2}i-1))}.$
(2)

Dividing the two equations leads to a standard constraint for oblique
pulsators with quadrupole modes:

$\tan i\tan\beta=4\frac{A_{+2}+A_{-2}}{A_{+1}+A_{-1}}.$ (3)

We can calculate the rotation inclination $i$ and magnetic obliquity $\beta$
of a quadrupole pulsator. Although this relation applies in the case of a pure
quadrupole mode, the results can provide us some information about the
geometry of the mode in HD 86181 for the pure case.

The determination of $\tan i\tan\beta$ for a dipole mode is similar to that
shown in eqn 3, but it is not possible to constrain $i$ and $\beta$
independently. However, for a normal quadrupole mode, eqns 1 and 2 provide two
equations in two unknowns, allowing us nearly uniquely to derive values for
$i$ and $\beta$. From eqns 1 and 3 we find $i=84\pm 3^{\circ}$, $\beta=30\pm
3^{\circ}$, or vice versa, for HD 86181. The uncertainties were calculated
through MCMC fits. With $i$, together with the rotation period and the
estimated radius, $v\sin i=6.5$ km s-1 can be derived. Although there is no
published $v\sin i$ value for comparison, this value is reasonable for a roAp
star.

For an axisymmetric quadrupole mode, the pulsation amplitude at the poles is
twice that at the equator and in antiphase. Maximum pulsation amplitude for
the angles determined above comes when $i-\beta=54^{\circ}$. Since the surface
nodes for an $\ell=2,m=0$ quadrupole lie at co-latitudes $\pm 54.7^{\circ}$,
at the time of pulsation maximum the pole is inclined $i-\beta=54^{\circ}$ to
the line of sight, one surface node is tangent to the lower limb of the star,
and the other surface node is over the the top limb. Hence we are seeing only
the pulsation polar cap at that time. Half a rotation later, the pole is
inclined by $i+\beta=114^{\circ}$; i.e., the pole we were seeing is now on the
other side of the star. The second pole has come into view, but is at poorer
viewing aspect, being inclined $66^{\circ}$ to the line of sight. That then
puts one of the surface nodes close to the line of sight, i.e.
$66-54=12^{\circ}$. Hence much of the visible hemisphere is dominated by the
equatorial region. Figure 7 shows schematically this geometry at four rotation
phases.

Figure 7: The schematic diagram of the viewing geometry of quadrupole mode of
HD 86181 through one rotation cycle. Red dots indicate the pulsation poles,
and blue dashed lines indicate the surface nodes at co-latitudes $\pm
54.7^{\circ}$.

For a pure oblique quadrupole mode, the pulsation amplitude distribution on
the surface, $A_{\theta}$, is proportional to $\frac{1}{2}(3\cos^{2}\theta-1)$
where $\theta$ is co-latitude, the angle to the poles. With knowledge of the
rotational inclination, $i$, and magnetic obliquity, $\beta$, we can calculate
an integral to obtain the pulsation amplitude at any time during a rotation
cycle. Numerically, the sphere surface of the star is divided into a grid;
then, with the formula, the pulsation amplitude for each cell of the grid can
be calculated. With $i$, $\beta$ and the rotation angle at a given time (t),
$2{\pi}{\nu}_{rot}t$, we know which grid cells can be seen by us and also the
projection of each cell. Then the integrated and projected surface pulsation
amplitude can be derived. The limb-darkening model for TESS (Claret, 2018) is
used here. The results are shown as the blue curves in Fig. 4. The maximum of
the integral pulsation amplitude is fixed to be the same as the one derived
from the model (red line). Since the calculation just considers the pulsation
as a pure quadrupole mode, the difference between blue and red line shows the
contribution from the radial and dipole component.

##  5 Spherical harmonic decomposition

Using the technique of Kurtz (1992), the quintuplet for HD 86181 can be
decomposed into a spherical harmonic series. This model is also based on the
oblique pulsator model. Although there are some caveats of this model,
estimates of pulsation amplitudes at some special phases and pulsation
amplitude ratios can be made easily.

The decomposition was done using the frequencies, amplitudes and phases from
Table 3. In order to interpret the two maximum pulsation amplitudes, we
calculated the decomposition with the time zero point $t_{0}={\rm
BJD}~{}2458569.26128$. The results are shown in Table 4.

Table 4: Results of the spherical harmonic decomposition (with the time zero point $t_{0}={\rm BJD}~{}2458569.26128$) of the quadrupole mode in HD 86181 for $i=84^{\circ}$ and $\beta=30^{\circ}$.  $\ell$ |  $A^{(\ell)}_{-2}$ (mmag) |  $A^{(\ell)}_{-1}$ (mmag) |  $A^{(\ell)}_{0}$ (mmag) |  $A^{(\ell)}_{+1}$ (mmag) |  $A^{(\ell)}_{+2}$ (mmag) |  $\phi$ (rad)
---|---|---|---|---|---|---
2 |  0.093 |  0.060 |  $-0.278$ |  0.056 |  0.081 |  -0.322
1 |  |  0.015 |  0.005 |  0.014 |  |  1.962
0 |  |  |  0.542 |  |  |  -0.204

In recent works, we have corrected a small error in the decomposition code.
The original code used to calculate the decomposition of HD 6532 (Kurtz et
al., 1996b) and several stars miscoded equations (8) and (10) in Kurtz (1992).

The decomposition components of HD 86181 show that at phase = 0, the dipole
$\ell=1$ component contributes only 0.034 mmag to the quadrupole mode – almost
nothing comparing to the strong radial contribution, which means that the
polar amplitude is increased and the equatorial amplitude is reduced compared
to a pure quadrupole mode. These results verify the assumption that the
pulsation amplitude maximum comes from the poles, with the secondary maximum
from the equator.

As an example, we estimate the pulsation amplitude maximum at phase 0.
According to the eqns (20), (21) and (22) in Kurtz (1992), at phase 0, the
pulsation amplitude is

$A=\sqrt{(\sum\limits_{\ell=0}^{2}\sum\limits_{m=-\ell}^{\ell}A^{\ell}_{m}\cos{\phi^{\ell}})^{2}+(\sum\limits_{\ell=0}^{2}\sum\limits_{m=-\ell}^{\ell}A^{\ell}_{m}\sin{\phi^{\ell}})^{2}}$.

The amplitude of the dipole mode ($A^{\ell=1}_{m}$) is negligible, and the
quadrupole and radial components have similar phases, meaning $\phi^{\ell=2}$
and $\phi^{\ell=0}$ can be considered as the same, so they can add at the time
of amplitude maximum. Therefore, the pulsation amplitude is
$A=0.542+0.093+0.060-0.278+0.056+0.081=0.559$ mmag, which fits the pulsation
phase plot well. Of course, the decomposition technique was designed to fit
the data, so it is not a surprise that it does. This discussion is to give a
mental picture of why this is so. More precisely, a fit of all three spherical
harmonic components taking into account that the exact phases seen in Table 4
gives the fit shown in Fig. 4 as the red curves.

In addition to the quintuplet for HD 86181, there are two doublets. With the
$i$ and $\beta$ in section 7, we derive $\tan i\tan\beta=4.76$. For dipoles
that gives

$\frac{A_{+1}+A_{-1}}{A_{0}}=4.76.$ (4)

We therefore expect to see triplets with very small central components at
$\nu_{2}$ and $\nu_{3}$, with amplitudes only about 0.03 mmag, which is at the
detection limit for these data. This supports the identification of $\nu_{2}$
and $\nu_{3}$ as dipole modes, and it is therefore no surprise that we see
doublets separated by twice the rotation frequency.

##  6 The large separation and acoustic cut-off frequency

The large separation, $\Delta\nu$, is the separation in frequency of modes of
the same degree and consecutive radial orders, and is proportional to the
square-root of the mean density of the star, i.e.,
$\Delta\nu\propto\sqrt{\rho}$ (e.g. Gabriel et al. 1985). This relation was
developed for the frequencies of high-order, acoustic, adiabatic, non-radial
oscillations (Tassoul, 1980, 1990). Since the roAp pulsations are in the
asymptotic regime, they are also applicable here. If the pulsation modes, or
at least relative radial orders are identified, the large separation can, in
principle, be determined.

To calculate the large separation, stellar radius and mass are required. We
estimate the radius of HD 86181 from $L=4{\pi}{\sigma}R^{2}T_{\rm eff}^{4}$,
and its mass from $M/{\rm M}_{\odot}=(L/{\rm L}_{\odot})^{1/4}$ (derived from
stellar homology relations; see, e.g., Eddington 1924) with the luminosity in
Table 1. We find $R=1.65$ R⊙, $M=1.72$ M⊙, and $\log g=4.19$ (cgs) for HD
86181. Although the mass is obtained from a rough scaling relation, it is fine
for estimating the large separation and the cut-off frequency in this section.

With the knowledge that the doublets we see are the result of dipole modes
with undetected central peaks, we are able to derive the mode frequencies to
be $\nu_{2}=230.103$ d-1 and $\nu_{3}=235.737$ d-1 by taking the average of
the two sets of sidelobes. That then gives the mode frequency separations to
be $\nu_{1}-\nu_{2}=2.668$ d-1 = 30.87 $\mu$Hz, and $\nu_{3}-\nu_{1}=2.967$
d-1 = 34.32 $\mu$Hz. Using the radius, mass and $\log$g estimated above and
the value of the solar large frequency separation
$\Delta\nu_{\odot}=134.88\,\mu$Hz (Huber et al., 2011), through
$\Delta\nu\propto\sqrt{\frac{g}{R}}$, we estimate $\Delta\nu/2=38.2$ $\mu$Hz,
which is consistent with the observations.

In roAp stars part of the pulsation mode energy can be refracted back into the
star by the influence of the magnetic field, even when the frequency of the
mode is above the acoustic cut-off frequency, $\nu_{ac}$ (Sousa & Cunha, 2008;
Quitral-Manosalva, Cunha & Kochukhov, 2018b). Therefore, there is no reason to
assume that very high frequency modes will not be observed in these pulsators.
Nevertheless, theory predicts that the excitation by the opacity mechanism
takes place in a frequency range that is close to, but does not exceed the
cut-off frequency and, thus, that an alternative excitation mechanism would be
required to excite modes of yet higher frequencies (Cunha et al., 2013). It is
therefore of interest to estimate the cut-off frequency in HD 86181 based on
the star’s global properties. Using the mass, radius and the effective
temperature in solar values in Table 1, and the scaling relation
$\nu_{ac}\propto g/\sqrt{T_{\rm eff}}$ (Brown et al., 1991) with
$\nu_{ac,\odot}=5.55$ mHz (Fossat et al., 1992), we find that in HD 86181
$\nu_{ac}\approx 3.03$ mHz, which is slightly larger than the observed mode
frequencies, around 2.73 mHz.

##  7 Modelling oblique quadrupole pulsations distorted by dipole magnetic
fields

In this section, we present comparisons of the observed amplitude and phase
modulations of HD 86181 with a quadrupole pulsation calculated by the method
of Saio (2005) including the effect of a dipole magnetic field. We assume that
the pulsations in roAp stars are axisymmetric with the pulsation axis aligned
with the axis of the dipole magnetic field. The strength of the field is
denoted by $B_{\rm p}$, the magnetic field strength at the poles.

In the presence of a magnetic field, the pulsation frequency is modified only
slightly (see Fig. 8), while the eigenfunction is distorted considerably
because the magnetic effect generates $\ell=0,4,6,\ldots$ components of
spherical harmonics in addition to the main $\ell=2$ component. (We have
included twelve components; i.e, up to $\ell=22$.) The eigenfunction gives
pulsation amplitude and phase on each point on the surface as a function of
the angle from the magnetic (or pulsation) axis. The amplitude/phase
distribution can be converted to observational amplitude/phase modulation as a
function of rotation phase (see Saio & Gautschy 2004 for details) for a set of
$(\beta,i)$. The method of comparison is also discussed in Shi et al. (2020).

According to the estimated luminosity range, we selected some models on the
1.65, 1.68, and 1.70 M⊙ evolutionary tracks as indicated by triangles in the
HR diagram of Fig. 9, in which the initial composition $(X,Z)=(0.70,0.02)$ is
adopted, while the helium abundance is assumed to be depleted to 0.01 (mass
fraction) in the layers above the second helium ionisation zone (polar model
in Balmforth et al. 2001). For a stellar model, we find, firstly without
including a magnetic field, a quadrupole mode having a pulsation frequency
close to $\nu_{1}=232.77\,{\rm d}^{-1}$. Then, we re-calculate the quadrupole
mode by taking into account the effect of an assumed dipole magnetic field of
$B_{\rm p}$.

For each case, an appropriate set of $(\beta,i)$ is determined by fitting the
amplitude modulation of HD 86181. Then, the phase modulation is compared with
the observations. Generally, for most assumed values of $B_{\rm p}$, the
obliquity and inclination angles $(\beta,i)$ can be determined easily by
fitting the predicted amplitude modulation with the observations, while the
theoretical phase modulation tends to be very small except for a certain range
of $B_{\rm p}$. Fig. 8 shows how theoretical phase modulations change with
changing $B_{\rm p}$ for a 1.68 M⊙ model. In this model, $6.5\lesssim B_{\rm
p}$/kG $\lesssim 8$ gives phase modulations that are comparable with the
observed ones. The required $B_{p}$ tends to be smaller in more massive models
because the mean density of the envelope is smaller in more massive stars.

Filled triangles in Fig. 9 indicate the loci of models whose amplitude and
phase modulations agree with the observed ones of HD 86181; agreements occur
if $B_{\rm p}\sim 9.0-6.0$ kG is assumed depending on the assumed stellar mass
of 1.65, 1.68, 1.70-M⊙. Among them, the three red triangles denote the models
whose large frequency separations agree with that of HD 86181. We have chosen
the 1.68-M⊙ model as the best model because the luminosity agrees with our
derived value better than the luminosity of the 1.70-M⊙ model does. However,
$\log T_{\rm eff}=3.859$ ($T_{\rm eff}=7230$ K) of the most approriate fit
model is somewhat lower than 7750 K listed in Table 1. This $T_{\rm eff}$
value is closer to $T_{\rm eff}=7320$ K obtained by McDonald, Zijlstra & Boyer
(2012) from a comparison of the SED with model atmospheres, and to $T_{\rm
eff}=7205$ K obtained by Masana, Jordi & Ribas (2006) from 2MASS photometry.

Fig. 10 compares amplitudes of the rotational sidelobes (top), amplitude
(middle) and phase (bottom) modulations between the best model with $B_{\rm
p}=7.0$ kG and HD 86181.By fitting the amplitude modulation, we find
$(\beta,i)$ = ($40^{\circ}$,$80^{\circ}$) each with an uncertainly of 5∘. The
$(\beta,i)$ given by the magnetically distorted model are only slightly
different from the pure quadrupole pulsator model: $i$ given by the distorted
model is consistent with the pure quadrupole pulsator model within the
1$\sigma$, while $\beta$ is consistent within 2$\sigma$. The range of the
phase modulation of the quadrupole model is small, which can be attributed to
contributions from $\ell=$4, 6, 8, …..

The dipole mode frequencies just above and below the quadrupole mode of the
best fitting model are 235.51 and 229.92 d-1, respectively, at $B_{\rm p}=7.0$
kG, which yield a large frequency spacing of 5.59 d-1 (or 64.7 $\mu$Hz), which
agrees with the observed large frequency spacing, $\nu_{3}-\nu_{2}=5.63$ d-1
(or 65.2 $\mu$Hz).222The frequency spacing of this model at $B_{\rm p}=0$ is
5.39 d-1 (or 62.4 $\mu$Hz).

For HD 42659, another roAp star pulsating in a distorted mode (Holdsworth,
Saio & Kurtz, 2019), the distorted model predicted the polar magnetic field
strength to be 0.8 kG by assuming that star pulsates in a quadrupole mode.
That result was consistent with the measured mean longitudinal magnetic field,
$\langle B_{l}\rangle=0.4$ kG (Kochukhov & Bagnulo, 2006; Hubrig et al.,
2006). However, the polar magnetic field strength predicted by our model for
HD 86181, $B_{p}=7.0$ kG, is significantly larger than the measured mean
longitudinal magnetic field, $\langle B_{l}\rangle=0.54$ kG (Bagnulo et al.,
2015). The cause of the difference is not clear. It could be a depth effect;
i.e., the magnetic field required in our model refers to the strength in the
hydrogen-rich envelope, while the measured magnetic field corresponds to the
strength in the outermost superficial layers. Also, there are some aspects
that are not considered in the model, such as the effects of surface spots.

Figure 8: Phase modulations (solid black lines) obtained by assuming various
strengths of magnetic fields for the quadrupole mode in the same model shown
in Fig. 10. Red dots are observed phase modulations of HD 86181, while dashed
magenta lines are the same as the one in the bottom panel of Fig. 10, which
are obtained from the oblique pulsator model of Kurtz (1992) (red lines in
Fig. 4). For all cases, $(\beta,i)=(40^{\circ},80^{\circ})$ are adopted, for
which the theoretical amplitude modulations are consistent with that of HD
86181, while Bp in the range 6.5 - 8 kG (e.g. 7 kG; Fig. 10) gives phase
modulation comparable with the observed one (see also Fig. 10). Figure 9: Loci
of roAp stars on the HR diagram with some evolutionary tracks with initial
composition of $(X,Z)=(0.70,0.02)$. The number along the ZAMS of each track
indicates the stellar mass in solar units. HD 86181 is shown in red and other
distorted quadrupole pulsators are shown in blue for comparison. (J1940 is not
shown because its location is very close to J1640.) Triangles on 1.70, 1.68
and 1.65 M⊙ tracks indicate the loci of models for which pulsation amplitude
and phase modulations are calculated; filled (open) triangles indicate models
whose phase modulations can (cannot) be fitted with the HD 86181 phase
modulation. Red filled triangles indicate models which have large frequency
spacings similar to the observed one. Parameters of roAp stars other than HD
86181 are adopted from Holdsworth et al. (2018b).  Figure 10: The amplitude
spectrum of rotational sidelobes (top panel) and amplitude (middle
panel)/phase (bottom panel) modulations of the quadrupole pulsation mode of HD
86181 are shown by red lines or dots. Dashed magenta lines (middle and bottom
panels) are obtained from the oblique pulsator model of Kurtz (1992) (red line
in Fig. 4). Black lines show the results of a best model of 1.68 M⊙ with
$B_{\rm p}=7$ kG, for which parameters are shown on the top of the diagram.

##  8 Driving of pulsations

The driving of pulsations in roAp stars is still a matter of debate. Non-
adiabatic pulsation calculations, assuming that envelope convection is
suppressed by the magnetic field at least in some angular region around the
magnetic pole, have been reasonably successful in explaining the driving of
most oscillations observed in roAp stars through the opacity mechanism acting
on the hydrogen ionization region (Balmforth et al., 2001; Cunha, 2002). The
same model also predicts that very high frequencies may be excited by the
turbulent pressure mechanism, a fact that has been suggested to explain the
pulsation frequencies observed in the roAp star $\alpha$ Cir (Cunha et al.,
2013). In this section we adopt the models discussed in these earlier works to
perform theoretical non-adiabatic pulsation calculations for HD 86181.

The analysis follows closely that presented by Cunha et al. (2013). In short,
the equilibrium model is derived from the matching of two spherically
symmetric models, one with envelope convection suppressed (the polar model)
and the other with convection treated according to a non-local mixing length
prescription (Spiegel, 1963; Gough, 1977a) (the equatorial model). It takes as
input the stellar mass, luminosity, effective temperature, chemical
composition (hydrogen, $X$, and helium, $Y$, mass fractions) and the
parameters associated with convection. The atmosphere is described by a
$T-\tau$ relation, which can be chosen amongst different options, with the
minimum optical depth, $\tau_{\rm min}$, being an additional input parameter.
Finally, helium settling can also be considered both in the polar and in the
equatorial regions, following a parametrized description with the surface
helium abundance in each region being additional input parameters.

The stability analysis is performed in each region separately and can consider
two different options for the surface boundary condition applied at the
minimum optical depth, namely, one that guarantees a full reflection of the
mode and one that allows waves with frequencies above the acoustic cut-off
frequency to propagate. In the equatorial model, the final non-adiabatic
solutions are computed using a non-local, time-dependent mixing-length
treatment of convection (Gough, 1977b; Balmforth, 1992). The results from the
non-adiabatic analysis in each region can then be combined to derive the
growth rates of modes in the model where convection is assumed to be
suppressed only in some angular region around the magnetic pole (the composite
model). Further details on the models can be found in Balmforth et al. (2001)
and references therein.

For each set of ($M$,$L$,$T_{\rm eff}$), four different physics configurations
were considered by varying different input parameters identified in previous
works as having significant impact on the stability results, namely: the
minimum optical depth, the outer boundary condition, and the amount of surface
helium. Table 5 summarizes the options in each case. Other parameters and
physics not mentioned here were fixed following the options adopted in
Balmforth et al. (2001).

Table 5: Modelling parameters for the cases illustrated in Fig. 11, all computed with $M=1.72\,{\rm M}_{\odot}$, $L=8.69\,{\rm L}_{\odot}$, $Y=0.278$, $X=0.705$. Model | Polar $Y_{\rm surf}$ | Equatorial $Y_{\rm surf}$ | $\tau_{\rm min}$ | Boundary | symbols
---|---|---|---|---|---
|  |  |  | condition | in Fig. 11
A | 0.01 | 0.278 | $3.5\times 10^{-5}$ | Reflective | circles
B | 0.01 | 0.278 | $3.5\times 10^{-4}$ | Reflective | squares
C | 0.01 | 0.278 | $3.5\times 10^{-5}$ | Transmissive | upward triangles
D | 0.1 | 0.1 | $3.5\times 10^{-5}$ | Reflective | rightward triangles

Fig. 11 shows an example of the results from the stability analysis in blue
and red, for polar and equatorial models, respectively, adopting the effective
temperature and luminosity in Table 1. Here we plot the relative growth rates
$\eta/\omega$ as a function of the cyclic pulsation frequency $\nu$, where
$\eta$ and $\omega$ are the imaginary and real parts of the angular
eigenfrequency, respectively, and a positive growth rate indicates the mode is
intrinsically unstable, thus excited. From the red symbols in the figure we
can see that all modes are stable in the equatorial model, independently of
the physics configuration adopted. In the polar models (blue symbols), a few
modes have positive growth rates at frequencies from $\sim$ 2.1 mHz up to
$\sim$ 2.7 mHz, depending on the physics considered. The range of excited
frequencies scales approximately with the square root of the mean density
(Cunha, 2002; Cunha et al., 2013). Given the uncertainty on the radius of the
star, one can thus confidently conclude that the region where the oscillation
frequencies are observed is within the range where the polar models predict
instability. Despite this, the growth rates on these polar models are one
order of magnitude smaller than the growth rates of the corresponding modes in
the equatorial model (in absolute value). This means that envelope convection
needs to be almost fully suppressed in order for these modes to be unstable in
the composite model (cf. figure 4 of Balmforth et al., 2001) and, thus,
explain the observations.

Figure 11: Normalized growth rates for polar (blue) and equatorial (red
symbols) regions as a function of the cyclic frequency $\nu=\omega/2\pi$.
Excited modes have positive growth rates. Different shape symbols represent
different modelling parameters we have used; circles are for model A in Table
5, while squares, upward triangles and rightward triangles for model B, C and
D, respectively. Zero growth rate is indicated by the horizontal dashed line
and the green shadowed region marks the range of observed frequencies.

##  9 Discussion and conclusions

We analysed HD 86181 with TESS data, and confirm it as a roAp star. The
rotation frequency is derived to be $\nu_{\rm rot}=0.48765\pm 0.00003$ d-1
($P_{\rm rot}=2.0507\pm 0.0001$ d). The pulsation frequency spectrum is rich,
consisting of one doublet, one quintuplet and another doublet. The central
frequency of the quintuplet is 232.7701 d-1 (2.694 mHz). The two doublets are
very likely to be sidelobes of two triplets, the amplitudes of whose central
frequencies are too small to be observed. With this interpretation, we
calculate the two central frequencies of the triplets to be 230.1028 d-1
(2.663 mHz) and 235.7361 d-1 (2.728 mHz).

Pulsation amplitude and phase modulation were calculated as a function of
rotation phase and shown to be modulated. Two maxima can be seen in the
rotational light curve, which indicates we see two primary spots in the TESS
pass-band, but the spot geometry is complex and further work is needed to
construct the chemical and magnetic map of this star.

We calculated the rotation inclination, $i$, and magnetic obliquity, $\beta$,
for HD 86181, which provided detailed information of the geometry and we used
those values with a spherical harmonic decomposition to better understand the
pulsation geometry and the distortion from a pure quadrupole mode.

Models considering the dipole magnetic field distortion were calculated and
compared with the observed amplitude and phase modulation. The best fit model
gives $B_{\rm p}=7.0$ kG and $(\beta,i)=(40^{\circ},80^{\circ})$. The
$(\beta,i)$ given by magnetic distortion model are only slightly different
from the pure quadrupole pulsator model with the relevant differences of (25
per cent, 5 per cent) for $(\beta,i)$, respectively. Also, the difference from
the phase modulation of the quadrupole model is small, which can be attributed
to higher degree components, $\ell=4,6,8,...$. The pulsation frequency and the
large frequency spacing given by this model are comparable with the
observation.

To explain the driving mechanism of this star, two non-adiabatic models were
constructed for HD 86181, one with envelope convection suppressed (the polar
model) and another considering convection (the equatorial model). We find that
polar model predicted the excitation of modes in the observed range.

The rich pulsation frequency spectrum let us study the large frequency
separation, $\Delta\nu$. The $\Delta\nu$ derived from $g$ and $R$ is
consistent with the observed value. The acoustic cut-off frequency,
$\nu_{ac}$, of this star is larger than the observed mode frequencies.

## References

  * Anders et al. (2019) Anders F. et al., 2019, A&A, 628, A94
  * Andrae et al. (2018) Andrae R. et al., 2018, A&A, 616, A8
  * Bagnulo et al. (2015) Bagnulo S., Fossati L., Landstreet J. D., Izzo C., 2015, A&A, 583, A115
  * Balmforth (1992) Balmforth N. J., 1992, MNRAS, 255, 603
  * Balmforth et al. (2001) Balmforth N. J., Cunha M. S., Dolez N., Gough D. O., Vauclair S., 2001, MNRAS, 323, 362
  * Balona, Holdsworth & Cunha (2019) Balona L. A., Holdsworth D. L., Cunha M. S., 2019, MNRAS, 487, 2117
  * Bigot & Dziembowski (2002) Bigot L., Dziembowski W. A., 2002, A&A, 391, 235
  * Bigot & Kurtz (2011) Bigot L., Kurtz D. W., 2011, A&A, 536, A73
  * Brown et al. (1991) Brown T. M., Gilliland R. L., Noyes R. W., Ramsey L. W., 1991, ApJ, 368, 599
  * Claret (2018) Claret A., 2018, A&A, 618, A20
  * Cunha (2002) Cunha M. S., 2002, MNRAS, 333, 47
  * Cunha et al. (2013) Cunha M. S., Alentiev D., Brandão I. M., Perraut K., 2013, MNRAS, 436, 1639
  * Cunha et al. (2019) Cunha M. S. et al., 2019, MNRAS, 487, 3523
  * Cunha, Fernandes & Monteiro (2003) Cunha M. S., Fernandes J. M. M. B., Monteiro M. J. P. F. G., 2003, MNRAS, 343, 831
  * Dziembowski & Goode (1985) Dziembowski W., Goode P. R., 1985, ApJ, 296, L27
  * Eddington (1924) Eddington A. S., 1924, Nature, 113, 786
  * Flower (1996) Flower P. J., 1996, ApJ, 469, 355
  * Fossat et al. (1992) Fossat E. et al., 1992, A&A, 266, 532
  * Freyhammer et al. (2009) Freyhammer L. M., Kurtz D. W., Elkin V. G., Mathys G., Savanov I., Zima W., Shibahashi H., Sekiguchi K., 2009, MNRAS, 396, 325
  * Gabriel et al. (1985) Gabriel M., Noels A., Scuflaire R., Mathys G., 1985, A&A, 143, 206
  * Gaia Collaboration (2020) Gaia Collaboration, 2020, VizieR Online Data Catalog, I/350
  * Gaia Collaboration et al. (2020) Gaia Collaboration, Brown A. G. A., Vallenari A., Prusti T., de Bruijne J. H. J., Babusiaux C., Biermann M., 2020, arXiv e-prints, arXiv:2012.01533
  * Gough (1977a) Gough D., 1977a, The current state of stellar mixing-length theory, Spiegel E. A., Zahn J. P., eds., Vol. 71, pp. 15–56
  * Gough (1977b) Gough D. O., 1977b, ApJ, 214, 196
  * Handler et al. (2006) Handler G. et al., 2006, MNRAS, 366, 257
  * Hauck & Mermilliod (1998) Hauck B., Mermilliod M., 1998, A&AS, 129, 431
  * Hey et al. (2019) Hey D. R. et al., 2019, MNRAS, 488, 18
  * Holdsworth (2021) Holdsworth D. L., 2021, Frontiers in Astronomy and Space Sciences, 8, 31
  * Holdsworth et al. (2021) Holdsworth D. L. et al., 2021, arXiv e-prints, arXiv:2105.13274
  * Holdsworth et al. (2018a) —, 2018a, MNRAS, 473, 91
  * Holdsworth et al. (2016) Holdsworth D. L., Kurtz D. W., Smalley B., Saio H., Handler G., Murphy S. J., Lehmann H., 2016, MNRAS, 462, 876
  * Holdsworth et al. (2018b) Holdsworth D. L., Saio H., Bowman D. M., Kurtz D. W., Sefako R. R., Joyce M., Lambert T., Smalley B., 2018b, MNRAS, 476, 601
  * Holdsworth, Saio & Kurtz (2019) Holdsworth D. L., Saio H., Kurtz D. W., 2019, MNRAS, 489, 4063
  * Holdsworth et al. (2018c) Holdsworth D. L., Saio H., Sefako R. R., Bowman D. M., 2018c, MNRAS, 480, 2405
  * Holdsworth et al. (2014) Holdsworth D. L., Smalley B., Kurtz D. W., Southworth J., Cunha M. S., Clubb K. I., 2014, MNRAS, 443, 2049
  * Huber et al. (2011) Huber D. et al., 2011, ApJ, 743, 143
  * Hubrig et al. (2006) Hubrig S., North P., Schöller M., Mathys G., 2006, Astronomische Nachrichten, 327, 289
  * Khomenko & Kochukhov (2009) Khomenko E., Kochukhov O., 2009, ApJ, 704, 1218
  * Kochukhov (2006) Kochukhov O., 2006, A&A, 446, 1051
  * Kochukhov (2009) —, 2009, Communications in Asteroseismology, 159, 61
  * Kochukhov & Bagnulo (2006) Kochukhov O., Bagnulo S., 2006, A&A, 450, 763
  * Kurtz (1982) Kurtz D. W., 1982, MNRAS, 200, 807
  * Kurtz (1985) —, 1985, MNRAS, 213, 773
  * Kurtz (1992) —, 1992, MNRAS, 259, 701
  * Kurtz et al. (1996a) Kurtz D. W., Marang F., van Wyk F., Roberts G., 1996a, MNRAS, 280, 1
  * Kurtz & Martinez (1994) Kurtz D. W., Martinez P., 1994, Information Bulletin on Variable Stars, 4013, 1
  * Kurtz et al. (1996b) Kurtz D. W., Martinez P., Koen C., Sullivan D. J., 1996b, MNRAS, 281, 883
  * Leckrone (1973) Leckrone D. S., 1973, ApJ, 185, 577
  * Masana, Jordi & Ribas (2006) Masana E., Jordi C., Ribas I., 2006, A&A, 450, 735
  * Mathys, Kharchenko & Hubrig (1996) Mathys G., Kharchenko N., Hubrig S., 1996, A&A, 311, 901
  * McDonald, Zijlstra & Boyer (2012) McDonald I., Zijlstra A. A., Boyer M. L., 2012, MNRAS, 427, 343
  * Montgomery & O′donoghue (1999) Montgomery M. H., O′donoghue D., 1999, Delta Scuti Star Newsletter, 13, 28
  * Moon & Dworetsky (1985) Moon T. T., Dworetsky M. M., 1985, MNRAS, 217, 305
  * Murphy, Shibahashi & Kurtz (2013) Murphy S. J., Shibahashi H., Kurtz D. W., 2013, MNRAS, 430, 2986
  * Napiwotzki, Schoenberner & Wenske (1993) Napiwotzki R., Schoenberner D., Wenske V., 1993, A&A, 268, 653
  * Perry (1991) Perry C. L., 1991, PASP, 103, 494
  * Pyper (1969) Pyper D. M., 1969, ApJS, 18, 347
  * Quitral-Manosalva, Cunha & Kochukhov (2018a) Quitral-Manosalva P., Cunha M. S., Kochukhov O., 2018a, MNRAS, 480, 1676
  * Quitral-Manosalva, Cunha & Kochukhov (2018b) —, 2018b, MNRAS, 480, 1676
  * Renson & Manfroid (2009) Renson P., Manfroid J., 2009, A&A, 498, 961
  * Saio (2005) Saio H., 2005, MNRAS, 360, 1022
  * Saio & Gautschy (2004) Saio H., Gautschy A., 2004, MNRAS, 350, 485
  * Shi et al. (2020) Shi F., Kurtz D., Saio H., Fu J., Zhang H., 2020, ApJ, 901, 15
  * Shibahashi & Takata (1993) Shibahashi H., Takata M., 1993, PASJ, 45, 617
  * Smalley et al. (2015) Smalley B. et al., 2015, MNRAS, 452, 3334
  * Sousa & Cunha (2011) Sousa J. C., Cunha M. S., 2011, MNRAS, 414, 2576
  * Sousa & Cunha (2008) Sousa S. G., Cunha M. S., 2008, MNRAS, 386, 531
  * Spiegel (1963) Spiegel E. A., 1963, ApJ, 138, 216
  * Stibbs (1950) Stibbs D. W. N., 1950, MNRAS, 110, 395
  * Takata & Shibahashi (1994) Takata M., Shibahashi H., 1994, PASJ, 46, 301
  * Takata & Shibahashi (1995) —, 1995, PASJ, 47, 219
  * Tassoul (1980) Tassoul M., 1980, ApJS, 43, 469
  * Tassoul (1990) —, 1990, ApJ, 358, 313
  * Trifonov et al. (2020) Trifonov T., Tal-Or L., Zechmeister M., Kaminski A., Zucker S., Mazeh T., 2020, A&A, 636, A74

## acknowledgements

This work was funded by the National Key R $\&$ D Program of China under grant
No.2019YFA0405504 and the National Natural Science Foundation of China (NSFC)
under grants No. 11973001, No. 11833002, No. 12090040 and No. 12090042. This
work includes data collected by the TESS mission. Funding for the TESS mission
is provided by the NASA Explorer Program. M. S. Cunha is supported by national
funds through FCT in the form of a work contract and through the research
grants UIDB/04434/2020, UIDP/04434/2020 and PTDC/FIS-AST/30389/2017, and by
FEDER - Fundo Europeu de Desenvolvimento Regional through COMPETE2020 -
Programa Operacional Competitividade e Internacionalização (grant:
POCI-01-0145-FEDER-030389). Daniel L. Holdsworth acknowledges financial
support from the Science and Technology Facilities Council (STFC) via grant
ST/M000877/1. Gerald Handler gratefully acknowledges funding through NCN grant
2015/18/A/ST9/00578. We thank the anonymous referee for a thorough,
knowledgeable review that improved this paper.

## Data Availability

The data underlying this article will be shared on reasonable request to the
corresponding author.