Tidal Effects on the Radial Velocities of V723 Mon: Additional Evidence for a Dark 3 M_⊙ Companion

We divide the stellar surface into 768 pixels ⁵ with equal solid angle ΔΩ using the HEALPix/healpy package (Górski et al. 2005; Zonca et al. 2019). ⁶ For each pixel labeled by j, we compute:

• 1.  
  normalized distance from the star's center R_j /R_⋆,
• 2.  
  normalized surface gravity g_j /g_⋆,
• 3.  
  foreshortening factor $\cos {\gamma }_{j}$, where γ_j is the angle between the surface normal and our line of sight,
• 4.  
  angle δ_j between the surface normal and radius vector,
• 5.  
  intensity ${I}_{j}({g}_{j},\cos {\gamma }_{j})$ that depends on g_j and $\cos {\gamma }_{j}$ through gravity- and limb-darkening, respectively,
• 6.  
  line-of-sight velocity with respect to the star's center of mass normalized by the equatorial rotation velocity, V_j /(2π R_⋆/P)

as a function of the orbital phase, orbital/spin inclination i, RG mass M_⋆, companion mass M_•, and semimajor axis scaled by the RG radius a/R_⋆. The flux contribution of each pixel ΔF_j is given by

$$\begin{eqnarray}&&{\rm{\Delta }}{F}_{j}\propto {I}_{j}({g}_{j},\cos {\gamma }_{j})\cos {\gamma }_{j}\displaystyle \frac{{R}_{j}^{2}{\rm{\Delta }}{\rm{\Omega }}}{\cos {\delta }_{j}}\end{eqnarray} \tag{ 1 }$$

for $\cos {\gamma }_{j}\gt 0$ (i.e., visible to the observer), and 0 otherwise. The flux change due to Doppler beaming is of order 10⁻⁴ for the rotation velocity of ≈20 km s⁻¹ and is not included. We also ignore the effects of irradiation and reflection because the companion appears to be nonluminous. A potential microlensing effect due to the compact companion eclipsing the RG star is also negligible given the large RG radius and relatively tight orbit (Trimble & Thorne 1969).

The quantities R_j , g_j , $\cos {\gamma }_{j}$, and $\cos {\delta }_{j}$ were computed assuming that the RG surface is described by a surface of the constant Roche potential (Wilson 1979), where R_j was solved iteratively for each grid point. The formulation is thus similar to the one in the PHOEBE model (e.g., Prša et al. 2016). The normalizations R_⋆ and g_⋆ were chosen to be the values at the points on the stellar equator perpendicular to the star-companion axis. The intensity I_j was computed adopting the quadratic limb-darkening law $I{(\cos \gamma )\propto 1-{u}_{1}(1-\cos \gamma )-{u}_{2}(1-\cos \gamma )}^{2}$, multiplied by ${({g}_{j}/{g}_{\star })}^{y}$ (Kopal 1959). Since the limb- and gravity-darkening coefficients u₁, u₂, and y are chromatic, the values of the coefficients need to be evaluated for an appropriate wavelength band, as will be discussed in Section 2.4. Note that the quadratic law is adopted considering a balance between the accuracy and computational cost, and that the model can be modified to incorporate more complex profiles as we try in Section 3.2. Although the quadratic law may fail to reproduce the limb-darkening at the very edge of the stellar disk accurately, the results from the modeling in this paper were found to be insensitive to the adopted profile.

One simple method to evaluate tidal RV anomalies is to compute the mean of V_j weighted by the flux ΔF_j , (∑_j V_j ΔF_j )/(∑_j ΔF_j ). This is the prescription proposed in the seminal works by Sterne (1941) and Wilson & Sofia (1976), and has been adopted in many other works. In reality, however, the RVs are derived via a more complicated procedure specific to each pipeline, and the simple "flux-weighted mean velocity" has been shown to deviate from actual measurements in the case of the Rossiter–McLaughlin signal (McLaughlin 1924; Rossiter 1924) originating from line-profile distortion due to transiting exoplanets (Ohta et al. 2005; Winn et al. 2005; Hirano et al. 2010, 2011). This is also found to be the case in our current problem: Figure 2 compares the tidal RV signal evaluated as the flux-weighted mean (crosses) against the values from a model described below (thick solid line), for the same set of model parameters (mean of the posterior distribution) from our analysis in Section 3. Because the RV data modeled in Section 3 were derived by computing cross-correlation between the observed spectra and a synthetic template spectrum (Strassmeier et al. 2012), we try to replicate the process as possible to model tidal RVs. The formulation here largely follows the one in Hirano et al. (2011).

We mainly consider the distortion of a single line at some specific wavelength and evaluate how this affects the RV values derived from a cross-correlation analysis. The line profile in velocity space ${ \mathcal F }(v)$, in the presence of rigid rotation and macroturbulence, is given by the following convolution

$$\begin{eqnarray}&&{ \mathcal F }(v)=S(v)\ast M(v).\end{eqnarray} \tag{ 2 }$$

Here

$$\begin{eqnarray}&&M(v)=\displaystyle \frac{{\displaystyle \sum }_{j}{{\rm{\Theta }}}_{j}(v-{V}_{j}){\rm{\Delta }}{F}_{j}}{{\displaystyle \sum }_{j}{\rm{\Delta }}{F}_{j}}\end{eqnarray} \tag{ 3 }$$

is the broadening kernel, where

$$\begin{eqnarray}\,\begin{array}{rcl}{{\rm{\Theta }}}_{j}(v) & = & \displaystyle \frac{1}{2}\left[{ \mathcal N }\,\left(v;0,\displaystyle \frac{1}{2}\,{\zeta }^{2}{\cos }^{2}{\gamma }_{j}\right)\right.\\ & & \left.+{ \mathcal N }\,\left(v;0,\displaystyle \frac{1}{2}\,{\zeta }^{2}{\sin }^{2}{\gamma }_{j}\right)\right],\\ { \mathcal N }(x;\mu ,{\sigma }^{2}) & \equiv & \displaystyle \frac{1}{\sqrt{2\pi {\sigma }^{2}}}\,\exp \left[-\displaystyle \frac{{\left(x-\mu \right)}^{2}}{2{\sigma }^{2}}\right]\end{array}\,\end{eqnarray} \tag{ 4 }$$

is the macroturbulence kernel for the radial-tangential model (Gray 2005). Here we assume equal contributions from the radial and tangential motions, and ignore the small Doppler shift due to flux difference between the rising and sinking gas streams (convective blueshift). ⁷ We assume that the intrinsic line profile S(v) in velocity space is given by a Gaussian ⁸

$$\begin{eqnarray}&&S(v)={ \mathcal N }(v;0,{\beta }_{{\rm{S}}}^{2}),\end{eqnarray} \tag{ 5 }$$

where ${\beta }_{{\rm{S}}}^{2}={\beta }_{\mathrm{thermal}}^{2}+{\beta }_{\mathrm{mic}}^{2}+{\beta }_{\mathrm{IP}}^{2}$ includes broadening contributions from the thermal motion, microturbulence, and instrumental profile, respectively. We assume β_thermal = 0.82 km s⁻¹ (corresponding to T_eff = 4500 K and iron atoms), β_mic = 1.0 km s⁻¹ (Holtzman et al. 2018; Jayasinghe et al. 2021), and β_IP = 2.31 km s⁻¹ (corresponding to the wavelength resolution R = 55,000) for V723 Mon. These yield the total β_S = 2.65 km s⁻¹. The resulting line profile is

$$\begin{eqnarray}\,\begin{array}{rcl}{ \mathcal F }(v) & \propto & \displaystyle \sum _{j}\left[{ \mathcal N }\,\left(v-{V}_{j};0,{\beta }_{{\rm{S}}}^{2}+\displaystyle \frac{1}{2}\,{\zeta }^{2}{\cos }^{2}{\gamma }_{j}\right)\right.\\ & & \left.+{ \mathcal N }\,\left(v-{V}_{j};0,{\beta }_{{\rm{S}}}^{2}+\displaystyle \frac{1}{2}\,{\zeta }^{2}{\sin }^{2}{\gamma }_{j}\right)\right]\,{\rm{\Delta }}{F}_{j}.\end{array}\,\end{eqnarray} \tag{ 6 }$$

The cross-correlation function (CCF) C(v) is computed by convolving a template spectrum T(v) against ${ \mathcal F }(v)$. We assume that the template is a theoretical spectrum similar to the observed one (as was the case in Strassmeier et al. 2012), but without broadening due to instrumental profile:

$$\begin{eqnarray}\,\begin{array}{rcl}T(v) & = & { \mathcal N }(v;0,{\beta }_{{\rm{T}}}^{2}),\\ {\beta }_{{\rm{T}}}^{2} & \equiv & {\beta }_{\mathrm{thermal}}^{2}+{\beta }_{\mathrm{mic}}^{2}.\end{array}\,\end{eqnarray} \tag{ 7 }$$

Then the CCF is given by

$$\begin{eqnarray}\,\begin{array}{rcl}C(v) & = & (T\star { \mathcal F })(v)\propto \displaystyle \sum _{j}\left[{ \mathcal N }\left(\,,v-{V}_{j};0,{\beta }_{{\rm{S}}}^{2}\right.\right.\\ & & \left.+{\beta }_{{\rm{T}}}^{2}+\displaystyle \frac{1}{2}\,{\zeta }^{2}{\cos }^{2}{\gamma }_{j}\right)+{ \mathcal N }\,\left(v-{V}_{j};0,{\beta }_{{\rm{S}}}^{2}\right.\\ & & +\left.\left.{\beta }_{{\rm{T}}}^{2}+\displaystyle \frac{1}{2}\,{\zeta }^{2}{\sin }^{2}{\gamma }_{j}\right)\right]\,{\rm{\Delta }}{F}_{j}.\end{array}\,\end{eqnarray} \tag{ 8 }$$

Thus the shape of the CCF profile and the resulting RVs do not only depend on V_j and ΔF_j , but on the line-profile parameters ζ and $\beta \equiv \sqrt{{\beta }_{{\rm{S}}}^{2}+{\beta }_{{\rm{T}}}^{2}}$.

The model RVs v_tidal are then derived as ${v}_{\mathrm{tidal}}\,={\mathrm{argmax}}_{v}C(v)$. Given that the derivatives of C(v) can be computed easily, v_tidal can often be found efficiently as the root of dC(v)/dv = 0 using the Newton–Raphson method. However, we found that this method sometimes fails when the tidal deformation is large and the CCF has many points of inflection. ⁹ Thus we adopted a slower but more robust procedure: we first cast dC(v)/dv = 0 into the form v = f(v) and solve it iteratively starting from the flux-weighted mean of V_j (which can be readily computed from the quantities in Section 2.1), and use two Newton–Raphson steps to make the iterative solution converge efficiently to the best one.

We have so far evaluated v_tidal by modeling the CCF for a single line at a particular wavelength, using the weights ΔF_j evaluated at single effective wavelength (see Section 2.1). In reality, however, the CCF is computed from the spectrum with many absorption lines at different wavelengths, and so the actual CCF would be a weighted sum of the many CCFs computed in Equation (8). Assuming that β and ζ are achromatic, this summation is equivalent to replacing ΔF_j in Equation (8) with the value integrated over the wavelength range of the spectrum with a certain weight W(λ). Since ΔF_j depends on the wavelength only through the limb- and gravity-darkening, this operation reduces to choosing the limb- and gravity-darkening coefficients u₁, u₂, and y evaluated in the appropriate band. The detailed shape of the weight W(λ) is difficult to model, because it depends on the strengths and amounts of the lines as well as the échelle orders used for RV extraction that affect the wavelength region to which the CCF, and hence the derived RV, is most sensitive. We thus introduce this effective band as another parameter that may vary within a physically reasonable range, and take into account its uncertainty in evaluating the coefficients. See Section 2.4 for practical implementation of this model.

We note that the formulation presented here is not necessarily a unique one but needs to be adjusted depending on the exact procedure adopted to extract RVs. For example, RVs may be derived from features of the CCF other than the peak or by fitting a Gaussian to the CCF; the CCF may be computed using a binary mask rather than a theoretical template spectrum; or the RVs may be derived by directly fitting the observed spectra with a theoretical model. Nevertheless the framework presented here remains useful for constructing similar models for RVs from different pipelines.

The RV measured at time t_i was modeled as

$$\begin{eqnarray}&&v({t}_{i})=-K\cos \left[2\pi \left(\displaystyle \frac{{t}_{i}-{t}_{0}}{P}\right)\right]+\gamma +{v}_{\mathrm{tidal}}({t}_{i}),\end{eqnarray} \tag{ 9 }$$

where t₀ is the time of the conjunction where the companion is in front of the red giant and γ denotes the RV zero-point. The RV semiamplitude K is

$$\begin{eqnarray}&&K={\left(\displaystyle \frac{2\pi G}{P}\right)}^{1/3}\,\displaystyle \frac{{M}_{\bullet }\sin i}{{\left({M}_{\star }+{M}_{\bullet }\right)}^{2/3}},\end{eqnarray} \tag{ 10 }$$

where G is Newton's gravitational constant.

We assume that the measurement errors for RVs are independent and identical Gaussians with the variance of ${\sigma }_{i}^{2}+{\sigma }_{\mathrm{jitter}}^{2}$. Here σ_i is an internal error of the ith data point, and σ_jitter models any other excess scatter that is not included in σ_i ¹⁰ and was inferred along with the other model parameters. Therefore the log-likelihood for a set of RV measurements y_i at times t_i is given by

$$\begin{eqnarray}&&\mathrm{ln}{ \mathcal L }=-\displaystyle \frac{1}{2}\displaystyle \sum _{i}\left\{\displaystyle \frac{{\left[{y}_{i}-v({t}_{i})\right]}^{2}}{{\sigma }_{i}^{2}+{\sigma }_{\mathrm{jitter}}^{2}}+\mathrm{ln}\left[2\pi ({\sigma }_{i}^{2}+{\sigma }_{\mathrm{jitter}}^{2})\right]\right\}.\end{eqnarray} \tag{ 11 }$$

The whole code was implemented using JAX (Bradbury et al. 2018) and NumPyro (Bingham et al. 2018; Phan et al. 2019). We assumed the priors as described in Table 1 and Section 2.4, and obtained posterior samples for the parameters using Hamiltonian Monte Carlo (Duane et al. 1987; Betancourt 2017). We sampled until the resulting chains had the split $\hat{R}\lt 1.01$ (Gelman et al. 2014) for all the parameters.

Table 1. System Parameters from Our RV Modeling

	Median &68.3% HPDI	90% HPDI	Prior
red giant mass M⋆ (M⊙)	${0.82}_{-0.14}^{+0.13}$	[0.62, 1.07]	${ \mathcal U }(0.5,3.0)$
red giant radius R⋆ (R⊙)	${24.25}_{-0.89}^{+0.88}$	[22.74, 25.67]	${ \mathcal N }(24.0,0.9,15)$
mass ratio M•/M⋆	${3.58}_{-0.54}^{+0.38}$	[2.85, 4.37]	${{ \mathcal U }}_{\mathrm{ln}}(\exp (0),\exp (3))$
companion mass M• (M⊙)	${2.95}_{-0.17}^{+0.17}$	[2.68, 3.23]	⋯
semimajor axis over red giant radius a/R⋆	${4.142}_{-0.093}^{+0.091}$	[4.001, 4.301]	⋯
RV semiamplitude K ( km s−1)	${65.268}_{-0.052}^{+0.061}$	[65.17, 65.36]	⋯
binary mass function (M⊙)	${1.7264}_{-0.0041}^{+0.0049}$	[1.7188, 1.7337]	⋯
time of conjunction t0 (BJD − 2450000)	${5575.0653}_{-0.0071}^{+0.0073}$	[5575.0533, 5575.0767]	${ \mathcal U }(5574.5954,5575.5954)$
orbital period P (days)	${59.9376}_{-0.0010}^{+0.0009}$	[59.9359, 59.9392]	${ \mathcal N }(60,1,58)$
cosine of orbital inclination $\cos i$	${0.12}_{-0.12}^{+0.06}$	[0.00, 0.28]	${ \mathcal U }(0,1)$
profile width β ( km s−1)	${3.93}_{-0.90}^{+0.42}$	[2.96, 4.99]	${ \mathcal N }(2.95,1,2.95)$
macroturbulence velocity ζ ( km s−1)	${5.52}_{-1.36}^{+0.97}$	[3.78, 7.48]	${ \mathcal N }(5.3,1,1)$
effective wavelength λeff (nm)	${626}_{-112}^{+100}$	[475.3, 802.5]	${ \mathcal N }(635.0,123.5,388)$
limb-darkening coefficient u1	${0.64}_{-0.20}^{+0.15}$	[0.37, 0.94]	⋯
limb-darkening coefficient u2	${0.17}_{-0.12}^{+0.13}$	[−0.04, 0.39]	⋯
gravity-darkening coefficient y	${0.46}_{-0.10}^{+0.09}$	[0.31, 0.64]	⋯
RV zero-point γ ( km s−1)	${1.885}_{-0.030}^{+0.031}$	[1.835, 1.937]	${ \mathcal U }(-10,10)$
RV jitter σjitter ( km s−1)	${0.168}_{-0.030}^{+0.026}$	[0.121, 0.214]	${{ \mathcal U }}_{\mathrm{ln}}(\exp (-5),\exp (0))$
projected rotation velocity $v\sin i$ ( km s−1)	${20.18}_{-0.88}^{+0.76}$	[18.80, 21.52]	⋯

Note. Values listed here report the medians and 68.3%/90% highest posterior density intervals (HPDIs) of the marginal posteriors, which were found to be unimodel for all the parameters. Priors—${ \mathcal N }(\mu ,{\sigma }^{2},l)$ means the normal distribution centered on μ and variance σ²; when l is specified the normal distribution is truncated at the lower limit l. ${ \mathcal U }(a,b)$ and ${{ \mathcal U }}_{\mathrm{ln}}(a,b)$ are the uniform- and log-uniform probability density functions between a and b, respectively. Dots indicate the parameters that were computed from the samples of the "fitted" parameters whose priors were explicitly specified.

Download table as:  ASCIITypeset image

We adopted the priors summarized in Table 1. They are uninformative unless otherwise specified below. We note that these priors are independent from the information derived from ellipsoidal variations or eclipses of the Balmer lines.

RG mass M_⋆ and radius R_⋆—We adopt a prior uniform in [0.5 M_⊙, 3.0 M_⊙] for the RG mass. For the RG radius, we assume a Gaussian prior ${ \mathcal N }({R}_{\star };24.0\,{R}_{\odot },0.9\,{R}_{\odot })$ based on the value derived from the SED, Gaia EDR3 distance (Gaia Collaboration et al. 2020), and correction for nonstellar flux using measured dilution of the absorption lines (Jayasinghe et al. 2021). Our definition of R_⋆ (Section 2.1) is not exactly the same as that in Jayasinghe et al. (2021), but the difference in the resulting solution is significantly smaller than the prior uncertainty and so can be ignored.

Macroturbulence ζ—The macroturbulence ζ, along with rotation, shapes the rotation kernel (Equation (4)) and can affect the RVs derived from CCFs. We adopt a Gaussian prior based on the relation from APOGEE DR13 (Holtzman et al. 2018), which encompasses the values estimated by Jayasinghe et al. (2021) from spectra and agrees with other measurements for RGB stars (e.g., Carney et al. 2008). The effect due to this uncertainty turns out to be minor, but would have been significant if the actual value was close to $v\sin i$ (see Section 3.1).

Profile width β—As detailed in Section 2.2, this parameter represents the broadening of the CCF due to intrinsic widths of the absorption lines and that of the template. Larger values of β tend to result in smaller v_tidal by smearing out the difference between pixels with different line-of-sight velocities V_j . The estimate in Section 2.2 gives β = 2.95 km s⁻¹, but the actual value could be larger depending on the factors including exact broadening of the theoretical template used for the analysis, microturbulence, and wavelength dependence of the instrumental profile. Thus we adopt a half-normal prior centered on 2.95 km s⁻¹ with the width of 1 km s⁻¹.

Limb- and gravity-darkening coefficients u₁ , u₂ , and y—Limb-darkening reduces the flux contribution from the surface elements with large line-of-sight velocities and reduces the amplitude of v_tidal. Gravity darkening, on the other hand, enhances the amplitude because the effect decreases the flux from the underrepresented side in velocity space (e.g., further reduces the "red" flux in Figure 1). Their values in our model depend on the effective wavelength band defined through the CCF (see Section 2.2). The effective band is difficult to quantify, but the following assumptions are reasonable: (i) the value is unlikely to be far from that obtained by integrating over the whole spectrum range with a uniform weight, because the lines relevant for RV measurements exist over the entire range, and (ii) the value should be bracketed by the values computed for the shortest and the longest wavelengths $({\lambda }_{\min },{\lambda }_{\max })$ in the spectrum. We implement this prior knowledge as follows. We take $u^{\prime} g^{\prime} r^{\prime} i^{\prime} z^{\prime}$-band coefficients theoretically computed with the ATLAS model (Claret & Bloemen 2011) for the effective temperature of 4500 K, log surface gravity of 1.5, and metallicity of −1 (Jayasinghe et al. 2021), and interpolate them over the central wavelengths of the bands to obtain the coefficients (u₁, u₂, y) as a function of wavelength. We then introduce a new parameter λ_eff that represents the effective band. This parameter was sampled from a Gaussian with the central value of $({\lambda }_{\min }+{\lambda }_{\max })/2$ and the width of $({\lambda }_{\max }-{\lambda }_{\max })/4$. Then we sample the coefficients from three independent Gaussians centered around u₁(λ_eff), u₂(λ_eff), and y(λ_eff) computed using the above deterministic relations and widths of 0.1 to incorporate uncertainties in the theoretical calculations. This prior is insensitive to the adopted spectroscopic parameters within their uncertainties as evaluated by Jayasinghe et al. (2021).

Projected rotation velocity $v\sin i$—Assuming tidal synchronization, our model automatically computes $v\sin i=2\pi {R}_{\star }\sin i/P$. This has also been evaluated from several different sets of spectra (see Jayasinghe et al. 2021), but we did not include this information in the fit for two reasons. First, interpretation of the "$v\sin i$" values of a tidally deformed star depends on how exactly they were extracted (Shahbaz 1998). Second, the measurements also depend on macroturbulence velocities adopted in those analyses, which are most likely different from each other and are not readily available. Nevertheless, the value predicted from our model turned out to be in reasonable agreement with those existing measurements.

Our RV model includes two additional parameters, macroturbulence ζ and profile width β, which are not required when the tidal RV is modeled as the flux-weighted mean velocity (Wilson & Sofia 1976). These parameters are not determined by the RV data but mostly determined by the adopted priors (see Table 1). Although we believe they are reasonable (and ζ is constrained from the spectra; Jayasinghe et al. 2021), we show in Figure 5 how the model v_tidal depends on these parameters to gauge their potential impacts on the other inferred parameters. Here the thick gray lines show the model computed for the mean parameter values of the posterior distribution, and the dashed and dotted lines show models where each parameter is perturbed by the values shown in the legends. The results show that it is essential to take into account their uncertainties as we did. In particular, the prior on the macroturbulence parameter needs to be chosen carefully. When the value is a significant fraction of $v\sin i$, as is the case for the dotted curve corresponding to ζ ≈ 11 km s⁻¹, the amplitude of the tidal RV depends significantly on this parameter. This can indeed be the case for some giant stars.

Zoom In Zoom Out Reset image size

The model atmospheres adopting spherical geometry suggest that the limb of giant stars with low surface gravity can be substantially darker than predicted by simple parametric laws (see, e.g., Figure 1 of Orosz & Hauschildt 2000) as we have adopted in the above analysis. To check on the sensitivity of our analysis to the limb-darkening profile, we redid the fit replacing the quadratic profile with the ones based on the intensity calculations from the PHOENIX model atmospheres with spherical geometry (Husser et al. 2013). ¹³ We used a model computed for the same atmospheric parameters as adopted in Section 2.4, retrieved the intensity values for 78 different $\cos \gamma$ and for the wavelengths spanning 352.5–952.5 nm at 100 nm intervals, and linearly interpolated them to compute $I(\cos \gamma ,{\lambda }_{\mathrm{eff}})$ in our RV model. We also incorporated 10% fractional uncertainty in the limb-darkening profile to take into account uncertainties in the model and the adopted atmospheric parameters. For this analysis, we adopted 3072 HEALPix pixels to ensure numerical stability for the updated limb-darkening profile with a sudden intensity drop at the limb. We found the results consistent with the above analysis using the quadratic law, including ${M}_{\bullet }={2.97}_{-0.17}^{+0.15}\,{M}_{\odot }$, $i={82.8}_{-3.3}^{+7.2}\,\deg$, and ${M}_{\star }={0.83}_{-0.15}^{+0.11}\,{M}_{\odot }$. Thus we conclude that our result is robust against the uncertainty of the limb-darkening profile. This appears reasonable given that the deviation from the quadratic law occurs mainly at $\cos \gamma \lesssim 0.25$, where the PHOENIX atmospheres predict almost zero intensities. The severe darkening results in the loss of stellar flux in the outermost ≈$1-\sqrt{1\mbox{--}{0.25}^{2}}=3 \%$ of the stellar disk. This is equivalent to a <1 km s⁻¹ change in $v\sin i$ and so plays a minor role in shaping the line profile. On the other hand, we found a slightly larger amplitude for the tidal RVs computed as the flux-weighted mean when the profile from the PHOENIX model was adopted.