Tidally driven mean flows in slowly and uniformly rotating massive main sequence stars

TL;DR

This study models tidally driven mean flows in rotating massive stars within binary systems, revealing how these flows depend on tidal frequency, resonance, and stellar parameters, with significant surface and equatorial effects.

Contribution

It introduces a linear oscillation-based method to compute tidally driven mean flows in rotating stars, highlighting the surface and equatorial localization of these flows.

Findings

Mean flows are dominated by the azimuthal component.

Flow amplitudes are significant near the surface where non-adiabatic effects matter.

Flow amplitudes depend on the mass ratio and orbital separation.

Abstract

We calculate tidally driven mean flows in a slowly and uniformly rotating massive main sequence star in a binary system. We treat the tidal potential due to the companion as a small perturbation to the primary star. We compute tidal responses of the primary as forced linear oscillations, as a function of the tidal forcing frequency $\omega_{\rm tide}=2(\Omega_{\rm orb}-\Omega)$, where $\Omega_{\rm orb}$ is the mean orbital angular velocity and $\Omega$ is the angular velocity of rotation of the primary star. The amplitude of the tidal responses is proportional to the parameter $f_0\propto (M_2/M)(a_{\rm orb}/R)^{-3}$, where $M$ and $M_2$ are the masses of the primary and companion stars, $R$ is the radius of the primary and $a_{\rm orb}$ is the mean orbital separation between the stars. For a given $f_0$, the amplitudes depend on $\omega_{\rm tide}$ and become large when $\omega_{\rm tide}$ is in resonance with natural frequencies of the star. Using the tidal responses, we calculate axisymmetric mean flows, assuming that the mean flows are non-oscillatory flows driven via non-linear effects of linear tidal responses. We find that the $\phi$-component of the mean flow velocity dominates. We also find that the amplitudes of the mean flows are large only in the surface layers where non-adiabatic effects are significant and that the amplitudes are confined to the equatorial regions of the star. Depending on $M_2/M$ and $a_{\rm orb}/R$, the amplitudes of mean flows at the surface become significant.

Tables

Table 1. Table 1: Complex eigenfrequency ω ¯ = ( ω ¯ R , ω ¯ I ) ¯ 𝜔 subscript ¯ 𝜔 R subscript ¯ 𝜔 I \bar{\omega}=(\bar{\omega}_{\rm R},\bar{\omega}_{\rm I}) of low radial order l = − m = 2 𝑙 𝑚 2 l=-m=2 g 𝑔 g -modes of the 15 ​ M ⊙ 15 subscript 𝑀 direct-product 15M_{\odot} ZAMS model for Ω ¯ = 0.1 ¯ Ω 0.1 \bar{\Omega}=0.1 .

	prograde		retrograde
mode	${\overline{\omega }}_R$	${\overline{\omega }}_I$	${\overline{\omega }}_R$	${\overline{\omega }}_I$
$g_1\cdots \cdots$	1.42585	2.39E-8	-1.45273	2.38E-8
$g_2\cdots \cdots$	0.91627	6.95E-8	-0.95505	8.24E-8
$g_3\cdots \cdots$	0.67962	2.36E-7	-0.72648	2.85E-7
$g_4\cdots \cdots$	0.53311	1.03E-6	-0.58549	1.00E-7
$g_5\cdots \cdots$	0.43741	5.66E-6	-0.49325	3.96E-6
$g_6\cdots \cdots$	0.37045	2.06E-5	-0.42877	1.63E-5
$g_7\cdots \cdots$	0.32112	8.12E-5	-0.38129	5.81E-5
$g_8\cdots \cdots$	0.28513	2.97E-4	-0.34638	2.21E-4
$g_9\cdots \cdots$	0.25756	6.64E-4	-0.32027	5.39E-4
$g_{10}\cdots \cdots$	0.23380	1.04E-3	-0.29800	3.56E-4

Equations

$${\partial\mbox{\boldmath$v$}^{\prime}\over\partial t}+2\mbox{\boldmath$\Omega$}\times\mbox{\boldmath$v$}^{\prime}=-{1\over\rho}\nabla p^{\prime}+{\rho^{\prime}\over\rho^{2}}\nabla p-\nabla\left(\Phi^{\prime}+\Psi\right),$$

$${\partial\rho^{\prime}\over\partial t}+\nabla\cdot\left(\rho\mbox{\boldmath$v$}^{\prime}\right)=0,$$

$${\nabla }^2{\Phi }^{\prime }=4\pi G{\rho }^{\prime },$$

$$\rho T\left({\partial s^{\prime}\over\partial t}+\mbox{\boldmath$v$}^{\prime}\cdot\nabla s\right)=\left(\rho\epsilon\right)^{\prime}-\nabla\cdot\mbox{\boldmath$F$}^{\prime},$$

$$\mbox{\boldmath$F$}^{\prime}_{\rm rad}=-\lambda_{\rm rad}^{\prime}\nabla T-\lambda_{\rm rad}\nabla T^{\prime},$$

$$\mbox{\boldmath$v$}^{\prime}={\rm i}\omega\mbox{\boldmath$\xi$},$$

$${\xi }_{r,\mathrm{e}}=-\frac{{\Phi }_{\mathrm{e}}^{\prime }+\Psi }{g},{\xi }_{h,\mathrm{e}}=\frac{1}{l(l+1)r}\frac{d}{dr}{r}^2{\xi }_{r,\mathrm{e}},$$

$${\rho }_{\mathrm{e}}^{\prime }=-{\xi }_{r,\mathrm{e}}\frac{d\rho }{dr},p_{\mathrm{e}}^{\prime }=-{\xi }_{r,\mathrm{e}}\frac{dp}{dr},$$

$${\nabla }^2{\Phi }_{\mathrm{e}}^{\prime }=4\pi G{\rho }_{\mathrm{e}}^{\prime },$$

$$\nabla\cdot\mbox{\boldmath$\xi$}_{\rm e}=0,$$

$${\nabla }^2\left({\Phi }_{\mathrm{e}}^{\prime }+\Psi \right)=\frac{4\pi G}{g}\frac{d\rho }{dr}\left({\Phi }_{\mathrm{e}}^{\prime }+\Psi \right),$$

$$\Psi =-f_0\frac{GM}{R}{x}^2{Y}_2^{-2}{e}^{2\mathrm{i}{\Omega }_{\mathrm{orb}}t}\equiv {\Psi }_2{Y}_2^{-2}{e}^{2\mathrm{i}{\Omega }_{\mathrm{orb}}t},$$

$$f_0=\sqrt{\frac{6\pi }{5}}\frac{G{M}_2/a_{\mathrm{orb}}^3}{{\sigma }_0^2}=\sqrt{\frac{6\pi }{5}}\frac{q}{{\left(a_{\mathrm{orb}}/R\right)}^3},$$

$$\frac{1}{x^2}\frac{d}{dx}\left(x^2\frac{d}{dx}F\right)-\left(\frac{6}{x^2}-k^2\right)F=0,$$

$$F=\frac{{\Phi }_{\mathrm{e}}^{\prime }+\Psi }{GM/R},k^2=-R^2\frac{4\pi G}{g}\frac{d\rho }{dr}.$$

$$F-\frac{x_0^2-k_0^2{x}_0^4/14}{2x_0-(2/7)k_0^2{x}_0^3}\frac{dF}{dx}=0,$$

$$\frac{d\ln \left|{\Phi }_e^{\prime }\right|}{d\ln r}=-3,$$

$$3F+\frac{dF}{dx}=-5f_0.$$

$$-\rho\omega^{2}\mbox{\boldmath$\xi$}+2{\rm i}\omega\Omega\rho\mbox{\boldmath$\Omega$}\times\mbox{\boldmath$\xi$}=-\nabla p^{\prime}-{\rho^{\prime}}\nabla\Phi-\rho\nabla\left(\Phi^{\prime}_{e}+\Psi\right).$$

$$\rho^{\prime}+\nabla\cdot\left(\rho\mbox{\boldmath$\xi$}\right)=0,$$

$${\rm i}\omega\rho T\delta s=\left(\rho\epsilon\right)^{\prime}-\nabla\cdot\mbox{\boldmath$F$}^{\prime},$$

$$\frac{{\rho }^{\prime }}{\rho }=-rA\frac{{\xi }_r}{r}+\frac{1}{{\Gamma }_1}\frac{p^{\prime }}{p}-{\alpha }_T\frac{\delta s}{c_p},$$

$$rA=\frac{d\ln \rho }{d\ln r}-\frac{1}{{\Gamma }_1}\frac{d\ln p}{d\ln r},$$

$${\Gamma }_1={\left(\frac{\partial \ln p}{\partial \ln \rho }\right)}_{\mathrm{ad}},{\alpha }_T=-{\left(\frac{\partial \ln \rho }{\partial \ln T}\right)}_p.$$

$${\xi }_r=r\sum _{j=1}^{j_{\max }}{S}_{l_j}(r)Y_{l_j}^m(\theta ,\varphi )e^{\mathrm{i}\omega t},$$

$${\xi }_{\theta }=r\sum _{j=1}^{j_{\max }}\left[H_{l_j}(r)\frac{\partial }{\partial \theta }{Y}_{l_j}^m(\theta ,\varphi )+T_{l_j^{\prime }}(r)\frac{1}{\sin \theta }\frac{\partial }{\partial \varphi }{Y}_{l_j^{\prime }}^m(\theta ,\varphi )\right]e^{\mathrm{i}\omega t},$$

$${\xi }_{\varphi }=r\sum _{j=1}^{j_{\max }}\left[H_{l_j}(r)\frac{1}{\sin \theta }\frac{\partial }{\partial \varphi }{Y}_{l_j}^m(\theta ,\varphi )-T_{l_j^{\prime }}(r)\frac{\partial }{\partial \theta }{Y}_{l_j^{\prime }}^m(\theta ,\varphi )\right]e^{\mathrm{i}\omega t},$$

$$p^{\prime }=\sum _{j=1}^{j_{\max }}{p}_{l_j}^{\prime }(r)Y_{l_j}^m\left(\theta ,\varphi \right)e^{\mathrm{i}\omega t},$$

$${\mathcal{T}}_2=-\int dV\overline{{\rho }_2^{\prime }\frac{\partial {\Phi }_2^{\prime }}{\partial \varphi }}=\frac{m}{2}{\int }_0^R\rho g{r}^{3}dr\mathrm{Im}\left(\frac{{\rho }_2^{\prime \ast }}{\rho }\frac{{\Phi }_{e,2}^{\prime }+{\Psi }_2}{gr}\right)\equiv \frac{G{M}^2}{R}{\overline{\mathcal{T}}}_2,$$

$$\overline{f}=\frac{1}{2\pi }{\int }_0^{2\pi }fd\varphi ,$$

Full text

††pubyear: 2015

1 Introduction

Tidal effects in binary systems of stars have long been investigated by many authors. The primary star in a binary system is affected by the gravitational field of the companion star that orbits around the primary, and vice versa. The tides affect binary evolution, leading to synchronization between the orbital motion and stellar rotation, circularization of the binary orbit, and change of the orbital separation between the stars (e.g., Hut, 1981). It was a common practice to consider equilibrium and dynamical tides separately. Equilibrium tides are tides considered in the limit of $\omega_{\rm tide}\rightarrow 0$, where $\omega_{\rm tide}$ is the forcing frequency caused by the orbital motion of the companion. Dynamical tides, on the other hand, are time dependent responses to the orbital motion of the companion. In the case of dynamical tides, frequency resonance between $\omega_{\rm tide}$ and natural frequencies of the star can take place and is expected to have significant effects on the binary evolution. It is dissipative processes accompanied by the tidal responses in binary stars that drive binary evolution.

Analytical and numerical studies of the tidal effects on binary evolution have been active since Zahn (1970, 1975, 1977) and Savonije & Papaloizou (1983, 1984) tried to estimate the time scales of synchronization and circularization of binary systems. In these studies, the tidal potential due to the companion star was assumed to be a small perturbation to the primary star so that the tidal responses of the primary are described by a linear theory of perturbations of stars. The magnitudes of the tidal responses is proportional to the parameter $f_{0}\propto(M_{2}/a_{\rm orb}^{3})(R^{3}/M)$, where $M$ and $R$ are the mass and radius of the primary star, $M_{2}$ is the mass of the companion star, and $a_{\rm orb}$ is the mean orbital separation between the stars. The tidal responses also depend on the forcing frequency $\omega_{\rm tide}$ and attain very large amplitudes when $\omega_{\rm tide}$ is in resonance with low frequency $g$-modes of the star. In a linear theory of perturbations, the resonant amplitudes of tidal responses are limited by dissipations such as produced by non-adiabatic effects and/or viscous effects accompanied with the responses. Although these early studies of tidal effects on binary systems of stars did not take account of the effects of stellar rotation on tidal responses, Savonije, Papaloizou, & Alberts (1995), Savonije & Papaloizou (1997), Witte & Savonije (1999ab, 2001, 2002), Ogilvie & Lin (2004, 2007) numerically investigated tidal responses of rotating stars. Stellar rotation brings about some complexities when estimating the tidal effects on binary stars and on the binary evolution. Because of the Coriolis force as a restoring force there appear rotational modes such as inertial modes and $r$-modes, whose frequencies are proportional to the rotation frequency $\Omega$ of the star (e.g., Unno et al 1989). Inertial modes propagate in isentropic regions found in the convective regions of stars, while $r$-modes, which are retrograde modes, propagate in the radiative envelope. If we consider tidal effects possibly caused by resonance between the forcing frequency $\omega_{\rm tide}$ and oscillation modes of rotating stars, we have to take into consideration rotational modes as well as $g$-modes when $\omega_{\rm tide}$ is comparable to or smaller than $\Omega$. Witte & Savonije (1999b, 2001), for example, discussed the effect of tidal locking on the binary evolution, and Ogilvie & Lin (2004) numerically investigated tidal excitation of inertial modes of a giant planet that has a large convective core and a thin radiative envelope.

Oscillations of rotating stars may excite axisymmetric mean flows in the stars. Lee et al. (2016) studied such mean flows driven by pulsationally unstable low frequency $g$- and $r$-modes of slowly pulsating B (SPB) stars, using a theory of wave-mean flow interaction (see Bühler 2014 for a review of the theory). In SPB stars, numerous low frequency oscillation modes are excited by the opacity bump mechanism operating in the temperature regions of $T\sim 1.5\times 10^{5}$K (e.g., Dziembowski et al 1993; Gautschy & Saio 1993). Lee et al. (2016) have shown that self-excited low frequency oscillations drive axisymmetric mean flows and that the $\phi$-component of the mean flow velocities dominates other components. The velocities of the mean flows become large in the surface layers of the envelope where non-adiabatic effects are significant. Note that, for mean flows driven by pulsationally unstable low frequency modes, the amplitudes of the modes and mean flows are undetermined within a linear theory of oscillation, unless we take account of amplitude saturation mechanisms such as non-linear couplings between oscillation modes (e.g., Lee 2012).

In this paper, we investigate axisymmetric mean flows driven by tidal responses of the primary star in a binary system. We treat the tidal responses, which are excited by orbital motion of the companion, as small amplitude perturbations of first order in the parameter $f_{0}$. We assume that axisymmetric mean flows of the second order are driven via non-linear effects of the responses. To compute mean flows for uniformly rotating stars, we use the formulation given by Lee et al. (2016), who employed an Eulerian perturbation theory of second order, where zonal averaging was used to pick up second order axisymmetric perturbations. We calculate tidal responses and mean flows as a function of the forcing frequency $\omega_{\rm tide}$.

We use a zero-age-main-sequence (ZAMS) star model of $15M_{\odot}$ as the background model for mean flow calculations. The ZAMS model has a chemically homogeneous and rather simple structure composed of a convective core and a radiative envelope and has a simple frequency spectrum of low frequency oscillation modes, which are all expected to be pulsationally stable. We also assume uniform and slow rotation of the star just for simplicity. For rapidly rotating stars, the frequency ranges of low radial order $g$-modes and inertial modes overlap, which would make the analyses more complicated. In §2, we give a brief description of the formulation we use for tidal response calculations and show some numerical results of the responses for the $15M_{\odot}$ ZAMS model. In §3, we describe numerical results for tidally driven axisymmetric mean flows. We conclude in §4.

2 Calculation of tidal responses

2.1 Basic equations for tidally perturbed stars

In a binary system of stars, the orbital motion of the companion star excites via tidal potential $\Psi$ time dependent tidal responses in the primary star, and vice versa. We let $\omega_{\rm tide}$ denote the forcing frequency associated with the tidal potential $\Psi$. If the tidal potential $\Psi$ is treated as a small perturbation to the primary star, the governing equations for tidal responses in the primary are given by a set of linearized basic equations of fluid dynamics:

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

where $\lambda_{\rm rad}={(4ac/3)}{(T^{3}/\kappa\rho)}$ is the radiative conduction coefficient, $v$ is the fluid velocity, $p$ is the pressure, $\rho$ is the mass density, $T$ is the temperature, $s$ is the specific entropy, $F$ is the energy flux, $\Phi$ is the gravitational potential, $\Psi$ is the tidal potential, $\epsilon$ is the nuclear energy generation rate per gram, $\kappa$ is the opacity, $G$ is the gravitational constant, $a$ is the radiation constant, $c$ is the velocity of light, and the primed quantities indicate Eulerian perturbations. Here, the companion star is assumed to be in the equatorial plane of the primary star, and the angular velocity of rotation $\Omega$ of the primary is assumed constant and parallel to the normal to the orbital plane. We have also assumed that the hydrostatic equilibrium in the primary star is given by $\nabla p=-\rho\nabla\Phi,$ that is, we have ignored equilibrium deformations due to rotation and tides. The energy flux $F$ is given by \mbox{\boldmathF}=\mbox{\boldmathF}_{\rm rad} in the radiative regions and \mbox{\boldmathF}=\mbox{\boldmathF}_{\rm rad}+\mbox{\boldmathF}_{\rm conv} in the convective regions, where \mbox{\boldmathF}_{\rm rad} and \mbox{\boldmathF}_{\rm conv} are the radiative and convective energy fluxes, respectively. For perturbations of the convective energy flux \mbox{\boldmathF}_{\rm conv}, we assume \delta\left(\nabla\cdot\mbox{\boldmathF}_{\rm conv}\right)=0 (see, e.g., Unno et al 1989), where $\delta$ indicates the Lagrangian perturbation.

In this paper, we assume that the time dependence of the perturbations is given by the factor $e^{{\rm i}\omega t}$ with $\omega$ being the oscillation frequency observed in the co-rotating frame of the star. For uniformly rotating stars, the Euler perturbations of the velocity, \mbox{\boldmathv}^{\prime}, is given by

[TABLE]

where \mbox{\boldmath\xi}=\xi_{r}\mbox{\boldmathe}_{r}+\xi_{\theta}\mbox{\boldmathe}_{\theta}+\xi_{\phi}\mbox{\boldmathe}_{\phi} is the displacement vector given in spherical polar coordinates $(r,\theta,\phi)$, and \mbox{\boldmathe}_{r}, \mbox{\boldmathe}_{\theta}, and \mbox{\boldmathe}_{\phi} are the orthonormal vectors in the $r$, $\theta$, and $\phi$ directions, respectively.

2.2 calculating equilibrium tide

As a response to the tidal potential $\Psi$, equilibrium tides may be defined as (e.g., Ogilvie & Lin 2004; see also Goldreich & Nicholson 1989)

[TABLE]

[TABLE]

and

[TABLE]

where $g=d\Phi/dr$. Note that, if we write \mbox{\boldmath\xi}_{\rm e}=\left(\xi_{r,\rm e}\mbox{\boldmathe}_{r}+\xi_{h,\rm e}\nabla\right)Y_{l}^{m}e^{{\rm i}\omega_{\rm tide}t}, we have

[TABLE]

indicating that the equilibrium tide is incompressible. Making use of equations (7), (8), and (9), we obtain

[TABLE]

where we have assumed $\nabla^{2}\Psi=0$. Integrating the differential equation (11) with appropriate boundary conditions, we obtain equilibrium tides $\Phi_{\rm e}^{\prime}$, \mbox{\boldmath\xi}_{\rm e}, $\rho_{\rm e}^{\prime}$ and $p_{\rm e}^{\prime}$.

We are interested in tidal responses excited by the potential given, in an inertial frame, as

[TABLE]

where $M$ and $R$ are the mass and radius of the primary star, $x=r/R$, $\Omega_{\rm orb}=\sqrt{G(M+M_{2})/a_{\rm orb}^{3}}$ is the mean angular velocity of the orbital motion, $a_{\rm orb}$ is the mean separation between the primary and companion stars, and

[TABLE]

where $M_{2}$ is the mass of the companion star, $q=M_{2}/M$, and $\sigma_{0}=\sqrt{GM/R^{3}}.$ For $l=2$, assuming $\Phi^{\prime}_{e}\propto Y_{2}^{-2}e^{2{\rm i}\Omega_{\rm orb}t}$, equation (11) may reduce to (see, e.g., Schwarzschild 1958)

[TABLE]

where

[TABLE]

We integrate the second order ordinary differential equation (14) from the centre to the surface of the star applying the boundary conditions given below. The boundary condition at the centre is the regularity condition given by

[TABLE]

where $x_{0}\ll 1$ and $k_{0}$ is the value of $k^{2}$ at $x=x_{0}$. The surface boundary condition at $x=1$ is given by (e.g., Ogilvie & Lin 2004)

[TABLE]

which leads to

[TABLE]

The amplitudes of the equilibrium and dynamical tides are proportional to the parameter $f_{0}$.

2.3 equations for tidal responses

Under the Cowling approximation, neglecting the Eulerian perturbation of the gravitational potential, we write the linearized equation of motion (1) as

[TABLE]

The perturbed continuity and entropy equations are

[TABLE]

[TABLE]

and the perturbed equation of state is given by

[TABLE]

where

[TABLE]

and

[TABLE]

Since separation of variables is not possible for the perturbations in rotating stars, we use finite series expansions in terms of spherical harmonic functions $Y_{l}^{m}(\theta,\phi)$ to represent the perturbations. Assuming that the equilibrium state is axisymmetric about the rotation axis, we expand the three components of the displacement vector \mbox{\boldmath\xi}(\mbox{\boldmathx},t) as

[TABLE]

[TABLE]

[TABLE]

and the Eulerian pressure perturbation, p^{\prime}(\mbox{\boldmathx},t), as

[TABLE]

where $l_{j}=2(j-1)+|m|$ and $l^{\prime}_{j}=l_{j}+1$ for even modes, and $l_{j}=2j-1+|m|$ and $l^{\prime}_{j}=l_{j}-1$ for odd modes for $j=1~{},2~{},3\cdots,~{}j_{\rm max}$ (e.g., Lee & Saio 1987). As indicated by the expressions given above, the perturbations are proportional to the factor $e^{{\rm i}m\phi+{\rm i}\omega t}$, and if we let $\sigma$ denote the oscillation frequency observed in an inertial frame, the oscillation frequencies $\omega$ in the co-rotating frame is given by $\omega=\sigma+m\Omega$. Substituting the series expansions of the perturbations into the perturbed basic equations (19) to (22) and (5), we obtain a finite set of linear ordinary differential equations for the expansion coefficients (see the Appendix). For a given tidal potential $\Psi$ and for a tidal forcing frequency $\omega=\omega_{\rm tide}$, we solve the finite set of differential equations with boundary conditions imposed at the centre and the surface of the star. The inner boundary conditions are the regularity condition for the perturbations and the condition for adiabatic oscillation given by $\delta s=0$. The outer boundary conditions are given by $\delta p=0$ and $\delta L_{\rm rad}=\delta\left(4\pi R^{2}\sigma_{\rm SB}T^{4}\right)$ with $\sigma_{\rm SB}$ being the Stefan-Boltzmann constant. See the Appendix for the detail.

2.4 calculating tidal responses

To investigate tidal responses of massive stars, we use a $15M_{\odot}$ zero age main sequence (ZAMS) model computed with a standard stellar evolution code using the OPAL opacity (Iglesias & Rogers 1996) for $X=0.7$ and $Z=0.02$. For this model, we plot $\psi_{2}/gr=(\Phi^{\prime}_{e,2}+\Psi_{2})/gr$ (solid line) and $\Psi_{2}/gr$ (dotted line) for $f_{0}=1$ in Figure 1. The difference between $\psi_{2}$ and $\Psi_{2}$ gives the contribution of the equilibrium tides $\Phi^{\prime}_{e}$, which become significant in the stellar core. Note that $(\xi_{r}/r)_{e}$ of the equilibrium tide is given by $(\xi_{r}/r)_{e}=-\psi_{2}/gr$.

We define the tidal torque ${\cal T}_{l}$ for $l=2$ as (e.g., Savonije & Papaloizou 1984)

[TABLE]

where

[TABLE]

and for a product of the first order perturbations $f_{1}$ and $f_{2}$, which are complex quantities, we may evaluate

[TABLE]

where the asterisk $(^{*})$ indicates complex conjugation.

In Figure 2 we plot the normalized torque $\left|\overline{{\cal T}_{2}}\right|$ as a function of the forcing frequency $\bar{\omega}_{\rm tide}$ for the $15M_{\odot}$ model for $\bar{\Omega}=0.1$ (left panel) and for $\bar{\Omega}=0.4$ (right panel), where $\bar{\omega}_{\rm tide}$ and $\bar{\Omega}$ denote dimensionless frequencies defined as $\bar{\omega}_{\rm tide}=\omega_{\rm tide}/\sigma_{0}$ and $\bar{\Omega}=\Omega/\sigma_{0}$. Since perturbations are assumed to be proportional to $e^{{\rm i}(\omega t+m\phi)}$ in this paper, the positive (negative) frequency $\omega_{\rm tide}$ corresponds to prograde (retrograde) forcing observed in the co-rotating frame of the star. As shown by the figure, there appear numerous peaks, produced by resonance between the forcing frequency $\omega_{\rm tide}$ and natural frequencies of $g$-modes and inertial modes of the star. The ZAMS model have a convective core and a radiative envelope, and $g$-modes propagate in the radiative envelope and inertial modes in the convective core where we have $N^{2}\approx 0$ with $N$ being the Brunt-Väisälä frequency. On the negative side of $\bar{\omega}_{\rm tide}$, we also find a sequence of resonance peaks associated with $r$-modes, which are retrograde modes propagating in the radiative envelope and have frequencies \bar{\omega}>\kern-11.99998pt\lower 4.73611pt\hbox{\sim}2m\bar{\Omega}/l^{\prime}(l^{\prime}+1)\approx-0.033 for $\bar{\Omega}=0.1$ and \bar{\omega}>\kern-11.99998pt\lower 4.73611pt\hbox{\sim}-0.1333 for $\bar{\Omega}=0.4$ when $m=-2$ and $l^{\prime}=3$. As shown by the figure, the tidal torque is significantly reduced in the inertial regime of $\left|\omega/\Omega\right|\leq 2$, except for that caused by the $r$-modes. For rapidly rotating stars, this inertial frequency regime overlaps the frequency ranges of low radial order $g$-modes. Although most of the conspicuous peaks result from resonance with $l=-m=2$ $g$-modes, we also find sequences of less pronounced peaks, which are produced by resonance with $g$-modes of $l=4$ and $m=-2$. For the $15M_{\odot}$ model, the tidal torques ${\cal T}_{2}$ has opposite signs between prograde and retrograde forcing and the sign stays the same as a function of $\omega_{\rm tide}$ except for very low frequency regions. For comparison, we tabulate the complex eigenfrequency $(\bar{\omega}_{\rm R},\bar{\omega}_{\rm I})$ of low radial order $l=-m=2$ $g$-modes of the $15M_{\odot}$ model for $\bar{\Omega}=0.1$, where $\bar{\omega}_{\rm I}>0$ indicates that the mode is pulsationally stable.

Figure 3 shows tidal responses of the $15M_{\odot}$ model at $\bar{\omega}_{\rm tide}=0.3704823794$ (left panel) and at $\bar{\omega}_{\rm tide}=0.35$ (right panel), where the real part of the expansion coefficient $S_{l}$ is plotted for $l=2$, 4, and 6. The left panel gives an example of tidal responses at a forcing frequency $\omega_{\rm tide}$ in resonance with a $l=-m=2$ $g$-mode, while the right panel shows a tidal response in off-resonance with low frequency modes, in which the response $S_{l=2}$ is approximately given by the equilibrium tide $(\xi_{r}/r)_{e}=-\psi_{2}/gr$ as shown by the long dashed line. In both cases, the component $S_{l=2}$ is dominating because the tidal potential $\Psi$ is here proportional to $Y_{l=2}^{m=-2}$. The amplitude at resonance can be much larger than the amplitude in off-resonance, which is comparable to $(\xi_{r}/r)_{e}$.

For comparison, we plot in Figure 4 the eigenfunction $S_{l}$ of the $l=-m=2$ $g_{6}$-mode and the derivative $dw/dr$ of the work function $w$ defined as (e.g., Unno et al 1989)

[TABLE]

where $\bar{\Omega}=0.1$ is assumed. We normalize the eigenfunction by $S_{l_{1}}=1$ at the surface Note that $dw/dr>0$ ($dw/dr<0$) indicates excitation (damping) regions for an oscillation mode. The $g_{6}$-mode is pulsationally stable, that is, the amount of damping exceeds that of driving in the interior. Comparing the left panels of Figures 3 and 4, we find that the tidal response $S_{l}$ at the resonance looks quite similar to the eigenfunction $S_{l}$, except for the amplitudes. The plot of $dw/dr$ in Figure 4 indicates that there extend an excitation region for the mode below $x=r/R\sim 0.94$, above which non-adiabatic damping prevails up to the surface.

3 Tidally Driven Mean Flows

We calculate mean flows driven by tidal responses, using the formulation given by Lee et al (2016). Here, tidal responses are considered as first order perturbations, while mean flows as second order perturbations in the parameter $f_{0}$.

3.1 perturbed equations of second order for mean flows

When tidal responses have small amplitudes and are regarded as a perturbation, any physical quantities f\left(\mbox{\boldmathx},t\right) of the primary star may be represented by

[TABLE]

where $f^{(0)}$ denotes the equilibrium quantities, $f^{(1)}$ the Euler perturbations of first-order, and $f^{(2)}$ the Eulerian perturbations of second-order in $f_{0}$. Similarly, the velocity field \mbox{\boldmathv}\left(\mbox{\boldmathx},t\right) may be expanded as

[TABLE]

and the equilibrium state is assumed to be that of a uniformly rotating star so that, in spherical polar coordinates $(r,\theta,\phi)$,

[TABLE]

where $\Omega$ is the angular velocity of rotation and assumed to be constant, and \mbox{\boldmathe}_{\phi} is the unit vector in the azimuthal direction. In this paper, we ignore equilibrium deformation caused by rotation and tidal force and assume that $f^{(0)}$ depends only on the radial distance $r$ from the center of the star. We apply the Cowling approximation to second order perturbations, neglecting the Euler perturbation $\Phi^{(2)}$ of the gravitational potential $\Phi$.

We employ a theory of wave-mean flow interaction to discuss axisymmetric flows driven by tidal responses in rotating stars (e.g., Bühler 2014). We regard the axisymmetric flows as mean flows, which contain both zero-th order and second order perturbations in $f_{0}$. The zero-th order quantities $f^{(0)}$ are those of equilibrium state, which is independent of time $t$. The first order quantities $f^{(1)}$ are tidal responses of the primary star. The second order quantities $f^{(2)}$ carry the time dependence of the mean flow. To derive governing equations for the second-order perturbations $f^{(2)}$ for mean flows, we use the zonal averaging defined by equation (30) and, assuming $\overline{f^{(1)}}=0,$ we obtain

[TABLE]

Here, we have ignored higher order terms $f^{(k)}$ with $k\geq 3$. Hereafter, we simply write $f^{(0)}$ and $f^{(2)}$ respectively for $\overline{f^{(0)}}$ and $\overline{f^{(2)}}$. The zonal averaging makes $f^{(0)}$ and $f^{(2)}$ independent of $\phi$.

We assume that non-oscillatory fluid flows arise in rotating stars via nonlinear effects of tidal responses $f^{(1)}$. Applying the zonal averaging to the basic equations, we obtain a set of differential equations that govern the second-order perturbations (Lee et al 2016):

[TABLE]

[TABLE]

[TABLE]

[TABLE]

where the first order quantities $f^{(1)}$ are simply written as $f^{\prime}$, and $\Phi^{(2)}$ and $\nabla^{2}\Phi^{(2)}=4\pi G\rho^{(2)}$ are ignored in the Cowling approximation.

For later convenience, we denote the right hand side of equation (37) as \mbox{\boldmathG}^{(2)}, that is,

[TABLE]

We may write \mbox{\boldmathG}^{(2)}, depending on the basis vector set, as

[TABLE]

where

[TABLE]

3.2 mean flow equations

The set of equations derived above are coupled linear partial differential equations for the second order axisymmetric perturbations where $(r,\theta)$ and $t$ are the independent variables. The products of first order perturbations provide inhomogeneous terms. To describe the $\theta$ dependence of the second order perturbations, we use series expansions of finite length, denoted by $k_{\rm max}$, in terms of spherical harmonic functions $Y_{l}^{0}(\theta,\phi)$. The velocity perturbation \mbox{\boldmathv}^{(2)} is expanded as

[TABLE]

[TABLE]

[TABLE]

and the pressure perturbation $p^{(2)}$ as

[TABLE]

where $l_{k}=2(k-1)$ and $l^{\prime}_{k}=l_{k}+1$ for $k=1,~{}2,~{}\cdots,~{}k_{\rm max}$.

By substituting the expansions (44) to (47) into equations (37) to (40), multiplying by a given spherical harmonic function, and integrating over solid angle, we derive a finite set of differential equations for the expansion coefficients, which depend on $r$ and $t$ (see Lee et al 2016). When we integrate over solid angle the non-linear terms such as (Y_{l_{k}}^{0})(\overline{\mbox{\boldmathv}^{\prime}\cdot\nabla\mbox{\boldmathv}^{\prime}})_{r}, we have to evaluate angular integration of products of three spherical harmonic functions, and we carry out the integration by introducing spin-weighted spherical harmonic functions ${}_{s}Y_{l}^{m}(\theta,\phi)$ (see, e.g., Newman & Penrose 1966; Varshalovich et al. 1988). Since the set of equations are linear equations for the second order expansion coefficients, we look for solutions whose time dependence is given by $e^{\gamma t}$. Replacing the time derivatives $\partial/\partial t$ by $\gamma$, the finite set of partial linear differential equations reduces to a set of linear ordinary differential equations that possess inhomogeneous terms.

Using vector notation, we formally write the set of linear ordinary differential equations with inhomogeneous terms as

[TABLE]

where

[TABLE]

and $\sf A$ and $I$ respectively represent the coefficient matrix and the inhomogeneous term (see Lee et al. 2016).

Non-radial components of equation (37) provide auxiliary equations, given by

[TABLE]

where

[TABLE]

with $\sf I$ being the identity matrix, and

[TABLE]

where $\bar{\gamma}=\gamma/\sigma_{0}$, $\bar{f}=\sqrt{4\pi/3}\bar{\Omega}$, and $\Lambda_{l}=l(l+1)$. See Lee et al (2016) for the definition of the matrices \mbox{\boldmath\sf\Lambda}_{0}^{1/2}, \mbox{\boldmath\sf C}_{B}^{0}, \mbox{\boldmath\sf C}_{B}^{1}, and \mbox{\boldmath\sf C}_{C}^{0}. The vectors \overline{\mbox{\boldmathG}_{q}^{j}} and \overline{\mbox{\boldmathG}_{\bar{q}}^{j}} for $j=0$ and $1$ on the right hand side of equation (50) come from the vector \mbox{\boldmathG}^{(2)} defined by equation (41). The $k$-th components of the vectors \overline{\mbox{\boldmathG}_{q}^{0}}, \overline{\mbox{\boldmathG}_{\bar{q}}^{0}}, \overline{\mbox{\boldmathG}_{q}^{1}}, and \overline{\mbox{\boldmathG}_{\bar{q}}^{1}} are given by

[TABLE]

where $do=\sin\theta d\theta d\phi$, and $l_{k}=2k-2$ and $l^{\prime}_{k}=l_{k}+1$ for $k=1,~{}2,~{}3,~{}\cdots$. We note that \mbox{\boldmath\sf W}=-\mbox{\boldmath\sf W}^{T} and $\sf W$ becomes singular when $\bar{\gamma}=0$. For non-zero values of $\bar{\gamma}$, we eliminate the variables \mbox{\boldmathz}_{h} and \mbox{\boldmathz}_{t} to derive equation (48) from the set of perturbed equations of second order, that is, we have to invert the matrix $\sf W$ for the elimination, which becomes numerically difficult when $|\bar{\gamma}|$ is extremely small.

The $\bar{\gamma}$ value may be determined by various processes. Dissipative effects such as non-adiabatic one in binary stars affect through tidal interactions the binary evolution which is described by slow changes of the binary parameters. The magnitude of the tidal effects on the primary is given by the parameter $f_{0}\propto qa_{\rm orb}^{-3}$, which depends on $a_{\rm orb}$ for a given value of $q$. The change rate of the mean separation $a_{\rm orb}$ may be given by (e.g., Savonije & Papaloizou 1997; Witte & Savonije 2002; Ogilvie & Lin 2004)

[TABLE]

where ${\cal T}_{2}$ is defined by equation (29), $E_{\rm orb}=-{GMM_{2}/2a_{\rm orb}}$, and we write the tidal forcing frequency as $\omega_{\rm tide}=n\Omega_{\rm orb}+m\Omega$ with $n=-m=2$. The normalized growth (decay) rate $\bar{\gamma}_{\rm tide}\sim 2d\ln f_{0}/d(\sigma_{0}t)=-6d\ln a_{\rm orb}/d(\sigma_{0}t)$ of the linear tidal responses may be given by

[TABLE]

where the sign of $\bar{\gamma}_{\rm tide}$ coincides with that of $\overline{{\cal T}_{2}}$. For $M_{2}\sim 0.1M$ and $a_{\rm orb}\sim 10R~{}{\rm to}~{}10^{2}R$, for example, we have $d\ln a_{\rm orb}/d(\sigma_{0}t)\sim-\overline{{\cal T}_{2}}$ and $f_{0}\sim 10^{-4}~{}{\rm to}~{}10^{-7}$. For this parameter range, even at resonance with $g$-modes, the amplitudes of the tidal responses $S_{l}$, which is proportional to $f_{0}$, is much smaller than unity (see Figure 3), and the magnitude of the normalized tidal torque $\overline{{\cal T}_{2}}$, which is proportional to $f_{0}^{2}$, is at most of order of $10^{-10}$ or smaller (see Figure 2). This suggests that $\bar{\gamma}_{\rm tide}$ is in general much smaller than $|2\bar{\omega}_{\rm I}|$ given in Table 1. If the orbital shrinkage takes place due to gravitational wave radiation, we may have (see, e.g., Landau & Lifshitz 1975)

[TABLE]

where $R_{g}=2GM/c^{2}$, and $\bar{\gamma}_{\rm grab}$ is generally smaller than $\bar{\gamma}_{\rm tide}$.

In this paper, we treat $\bar{\gamma}$ as a constant parameter of order of $10^{-8}$ so that we can properly inverse the matrix $\sf W$. We confirm that the flow patterns of tidally driven mean flows for $|\bar{\gamma}|=10^{-8}$ are the same as those for $|\bar{\gamma}|=10^{-10}$, and that for sufficiently small values of $\bar{\gamma}$, the magnitudes of $v_{\phi}^{(2)}$ scales as $\gamma v_{\phi}^{(2)}\approx{\rm constant}$ where the constant is almost independent of $\bar{\omega}_{\rm tide}$.

3.3 angular momentum transport by waves

The angular momentum transport by waves in rotating stars may be described by (e.g., Lee 2013, see also Grimshaw 1984)

[TABLE]

where $\xi$ is the displacement vector associated with the wave, and

[TABLE]

is the specific angular momentum in the $z$-direction, where

[TABLE]

[TABLE]

The total time derivative on the left-hand-side of equation (57) is defined as

[TABLE]

where

[TABLE]

and

[TABLE]

is the second order Lagrangian perturbation of the velocity, and the semicolon indicates covariant derivatives. Note that for uniform rotation, $v^{(0)}_{i;j;k}=0$ and hence \delta v_{i}^{(2)}(\mbox{\boldmathx})=v_{i}^{(2)}(\mbox{\boldmathx})+\overline{v_{i;j}^{(1)}(\mbox{\boldmathx})\xi_{j}(\mbox{\boldmathx})}. The treatment given above is based on the Lagrangian mean wave-mean flow interaction theory developed by Andrews & McIntyre (1976, 1978ab). See also Dunkerton (1980), Grimshaw (1984), and Bühler (2014) for reviews of wave-mean flow interaction theories.

Integrating equation (57) over the whole volume of the star, we obtain

[TABLE]

where dV=d^{3}\mbox{\boldmathx} and \tilde{\rho}(\mbox{\boldmathx})d^{3}\mbox{\boldmathx}=\rho(\hat{\mbox{\boldmathx}})d^{3}\hat{\mbox{\boldmathx}} with \hat{\mbox{\boldmathx}}=\mbox{\boldmathx}+\mbox{\boldmath\xi}(\mbox{\boldmathx}), and we have ignored the surface term assuming $p^{\prime}/\rho$ is finite at the surface, that is, $p^{\prime}/\rho$ remains finite as $\rho\rightarrow 0$ toward the surface. The right-hand-side becomes the tidal torque when we replace $\Phi^{\prime}$ by the tidal potential $\Phi^{\prime}_{e}+\Psi$.

If we assume $v_{r}^{(0)}=v_{\theta}^{(0)}=0$ for uniform rotation, equation (57) becomes

[TABLE]

where $\Phi^{\prime}=\Phi^{\prime}_{e}+\Psi\propto Y_{l=2}^{m=-2}$ for tidal responses. Integrating over solid angle, we obtain

[TABLE]

and it is convenient to denote the right hand side of equation (66) as ${\cal W}(r)$, that is,

[TABLE]

Equation (57) may be regarded as a mean flow equation that describes responses of the mean flow to waves. If the waves are non-dissipative, the right-hand-side of equation (57) vanishes, indicating conservation of the specific angular momentum \overline{\ell(\hat{\mbox{\boldmathx}})} as stated by Goldreich & Nicholson (1989). Responses of the waves to mean flows, on the other hand, may be described by the equation for wave action. In the Lagrangian mean theory (Andrwes & McIntyre 1978ab; see also Dunkerton 1980; Grimshaw 1984), the wave action $A$ obeys

[TABLE]

where

[TABLE]

and $K_{ij}$ is the $(i,j)$th cofactor of the Jacobian J\equiv{\rm det}(\partial\hat{\mbox{\boldmathx}}/\partial\mbox{\boldmathx}), and in small amplitude limit of the waves the dissipation term $D$ reduces to

[TABLE]

It is dissipative processes that cause interaction between the mean flows and waves.

3.4 velocity fields of tidally driven mean flows

For tidally driven mean flows, we calculate the velocity fields \mbox{\boldmathv}_{H}^{(2)}=v_{y}^{(2)}\mbox{\boldmathe}_{y}+v_{z}^{(2)}\mbox{\boldmathe}_{z} on three spherical surfaces of different radii $r/R=0.99$, 0.95 and 0.90, where, assuming the $x$-axis is towards the observer, the velocity fields $(v_{y},v_{z})$ in the $y$-$z$ plane are given by

[TABLE]

[TABLE]

and $\theta$ and $\phi$ are respectively the colatitude, measured from the $z$-axis, and the azimuthal angle, measured from the $x$-axis. For the mean flow calculations we use the expansion length $k_{\rm max}=16$, which is we find long enough to get good convergence of the expansions for the perturbations. Figure 5 shows \mbox{\boldmathv}_{H}^{(2)} at the forcing frequency $\bar{\omega}_{\rm tide}=0.3704823794$ that is in resonance with the prograde $l=-m=2$ $g_{6}$-mode of the $15M_{\odot}$ model for $\bar{\gamma}=-10^{-8}$. On each of the spherical surfaces, \mbox{\boldmathv}_{H}^{(2)} is normalized by its maximum value v_{\rm max}^{(2)}(x=r/R)\equiv{\rm max}(|\mbox{\boldmathv}_{H}^{(2)}(r,\theta,\phi)|) on that surface, and the length of the arrows is proportional to the magnitude of normalized \mbox{\boldmathv}_{H}^{(2)}. As discussed in Lee et al (2016), the $\phi$ component of \mbox{\boldmathv}_{H}^{(2)} dominates the $r$ and $\theta$ components, and hence the velocity fields \mbox{\boldmathv}_{H}^{(2)} of tidally driven mean flows are almost parallel to the equator of the star. The velocity field \mbox{\boldmathv}^{(2)}_{H} is symmetric about the equator, and the amplitudes of \mbox{\boldmathv}_{H}^{(2)} tend to be confined to the equatorial regions. Since mean flows arise from non-adiabatic effects accompanied with the responses, $v_{\rm max}^{(2)}(x)$ becomes largest in the outer most layers where the non-adiabatic effects become most significant. For example, we find $v_{\rm max}^{(2)}(x=0.9)\sim 0.01\times v_{\rm max}^{(2)}(x=0.99)$. In the equatorial regions, the velocities are prograde in the surface layers, while they become retrograde in the deep interior, suggesting that there arises differential rotation in radial direction. The amplitude confinement of \mbox{\boldmathv}_{H}^{(2)} into the equatorial regions also indicates differential rotation in the $\theta$-direction.

Figure 6 shows the velocity fields \mbox{\boldmathv}_{H}^{(2)} at the same forcing frequency $\omega_{\rm tide}$ but for $\bar{\gamma}=+10^{-8}$. The directions of the velocity \mbox{\boldmathv}_{H}^{(2)}\approx v_{y}^{(2)}\mbox{\boldmathe}_{y} are opposite to those for $\bar{\gamma}=-10^{-8}$, but we find that the directions of $v_{x}^{(2)}$ and $v_{z}^{(2)}$ remain the same. As indicated by equation (65), only the term $\partial\ell^{(2)}/\partial t=\gamma\ell^{(2)}$ explicitly depend on $\gamma$, and the terms other than $v_{\phi}^{(2)}$ in $\ell^{(2)}$ are products of the first order perturvations, which are assumed to be proportional to $e^{\gamma t/2}$. If $\partial v_{\phi}^{(2)}/\partial t=\gamma v_{\phi}^{(2)}$ is dominating, $v_{\phi}^{(2)}$ at a given $r$ has to change its sign according to the sign of $\gamma$ to balance the right-hand-side of equation (65), which does not explicitly depend on $\gamma$.

Figure 7 shows \mbox{\boldmathv}_{H}^{(2)} of mean flows at a forcing frequency $\omega_{\rm tide}$ in off-resonance with low frequency modes of the star, where we use $\bar{\omega}_{\rm tide}=0.35$ for $\bar{\gamma}=-10^{-8}$. The magnitudes of \mbox{\boldmathv}_{H}^{(2)} are much smaller than those at resonance with the $g_{6}$-mode. This is of course because the amplitudes of tidal responses in off-resonance are much smaller than those in the resonance. The velocity fields \mbox{\boldmathv}_{H}^{(2)} are confined to equatorial regions in the surface layers, but the confinement is not necessarily strong in the deep interior, where \mbox{\boldmathv}_{H}^{(2)} shows more complicated behavior as a function of $\theta$. The velocities \mbox{\boldmathv}_{H}^{(2)} are retrograde in the surface layers, but \mbox{\boldmathv}_{H}^{(2)} at the equator can be prograde at the surface. It may be important to note that the averaged velocities \int\sin\theta\mbox{\boldmathv}_{H}^{(2)}do in the deep interior are prograde, which may be consistent with the belief that prograde tidal forcing causes acceleration of rotation rate of the star.

Figure 8 plots the function ${\cal W}(x)$ (solid line) at the two tidal forcing frequencies $\bar{\omega}_{\rm tide}=0.3704823794$ and $\bar{\omega}_{\rm tide}=0.35$, where the dotted and dashed lines represent the first and second terms on the right-hand-side of equation (67), respectively. The two terms cancel each other to lead to small amplitude ${\cal W}$ in the deep interior. The function ${\cal W}$ has large amplitudes only in the outer layers of the envelope. The amplitudes $|{\cal W}|$ at the resonance is much larger than those in off-resonance. The $r$-dependence of $\cal W$ for the response in the $g$-mode resonance looks quite similar to that of the function $(\rho/\rho_{m})^{-1}dw/dr$ for the eigen $g$-mode. This similarity may suggest that in the case of resonant forcing the velocity fields $v_{\phi}^{(2)}$ of mean flows are closely related to the damping and driving regions for the oscillation mode. It is interesting to note that the function $\cal W$ for the off-resonance forcing behaves quite differently from that for the resonant forcing. The function $\cal W$ for off-resonance forcing is positive in the surface layers, while it is negative for the resonant forcing. If the term $\int do~{}\partial\ell^{(2)}/\partial t$ is dominating on the left-hand-side of equation (66) and the approximation $\int do~{}\partial\ell^{(2)}/\partial t\approx\int do~{}r\sin\theta\partial v_{\phi}^{(2)}/\partial t$ is valid, $\int do~{}\sin\theta\partial v_{\phi}^{(2)}/\partial t=\int do~{}\sin\theta\gamma v_{\phi}^{(2)}$ is positive (negative) where ${\cal W}$ is positive (negative), which is what we find for the prograde forcing.

Assuming $\bar{\gamma}=-10^{-8}$, we calculate the mean flow velocity \mbox{\boldmathv}_{H}^{(2)} for the retrograde forcing $\omega_{\rm tide}$ in resonance with the $l=-m=2$ $g_{7}$-mode (Figure 9) and in off-resonance with $g$-modes (Figure 10). The amplitudes of \mbox{\boldmathv}_{H}^{(2)} tend to be confined to the equatorial regions in the surface layers although this confinement becomes weaker in the deep interior, particularly for the off-resonant forcing. For the resonant forcing, the velocities \mbox{\boldmathv}_{H}^{(2)} are retrograde at the surface but become prograde as we go into the deep interior, which is consistent with the behavior of the function ${\cal W}(r)$ in the left panel of Figure 11. Note that the signs of the function $\cal W$ at a given radial distance $r$ are in general opposite to each other between the prograde and retrograde forcing with similar $|\bar{\omega}_{\rm tide}|$. For the off-resonant forcing, the velocities \mbox{\boldmathv}_{H}^{(2)} are mostly retrograde, which, as shown in the right panel of Figure 11, is not consistent with the interpretation in terms of ${\cal W}(r)$ based on the assumption that $\int do~{}\sin\theta\gamma v_{\phi}^{(2)}$ is dominating. This may suggest that the term $\int do~{}\sin\theta\partial v_{\phi}^{(2)}/\partial t$ is not necessarily dominating on the left hand side of equation (66) for off-resonance tidal forcing.

Figure 12 plots $v_{\rm max}^{(2)}/(R\sigma_{0})$ at $x=0.99$ as a function of the forcing frequency $\omega_{\rm tide}$ for the $15M_{\odot}$ model for $\bar{\gamma}=-10^{-8}$ and $f_{0}=1$. The velocity $v_{\rm max}^{(2)}/(R\sigma_{0})$ makes peaks at resonance with low frequency modes and can be as large as $10^{12}$ for low radial order $g$-modes, and the height of the peaks decreases as the radial order of $g$-modes increases. We also note that $v_{\rm max}^{(2)}/(R\sigma_{0})$ in off-resonance stays around $\sim 10^{6}$. Since $v_{\rm max}^{(2)}/(R\sigma_{0})\propto f_{0}^{2}\propto q^{2}(a_{\rm orb}/R)^{-6}$, if we assume $q\sim 0.1$, $v_{\rm max}^{(2)}/(R\sigma_{0})$ at $\bar{\omega}_{\rm tide}$ in resonance with low radial order $g$-modes will be $\sim 10^{4}$ for $a_{\rm orb}/R\sim 10$ and $\sim 10^{-2}$ for $a_{\rm orb}/R\sim 10^{2}$.

4 conclusion

In this paper, we computed tidally driven axisymmetric mean flows in a slowly and uniformly rotating massive main sequence star in a binary system, assuming that the tidal potential due to the companion star is a small perturbation to the primary star and the mean flows excited in the primary are of second order of the perturbation amplitudes. Here, we ignored equilibrium structure deformation caused by rotation and tidal force so that the equilibrium structure can be treated as being spherical symmetric. To compute the mean flows, we made a simplifying assumption that the time derivatives \partial\mbox{\boldmathv}^{(2)}/\partial t can be replaced by \gamma\mbox{\boldmathv}^{(2)} where $\gamma$ is a constant parameter regarded as the growth (or decay) rate of the second order perturbations. We find that the $\phi$-component of the velocity fields \mbox{\boldmathv}^{(2)} is the dominant one and that the amplitudes tend to be confined in the equatorial regions in the surface layers and decrease as we go into the deep interior where non-adiabatic effects become insignificant. We find the velocities \mbox{\boldmathv}^{(2)} in the deep interior are prograde (retrograde) for the prograde (retrograde) forcing $\omega_{\rm tide}$, which may be consistent with the picture that dissipation in the deep interior associated with tidal responses cause synchronization between orbital motion and stellar rotation in binary systems. We also discussed the relation between the term $\partial v_{\phi}^{(2)}/\partial t$ averaged over the colatitude $\theta$ and the function ${\cal W}(r)$, assuming the averaged $\partial v_{\phi}^{(2)}/\partial t$ is the dominant term in the angular momentum conservation equation.

The velocities $v_{\phi}^{(2)}$ of tidally driven mean flows depend on both $r$ and $\theta$, which inevitably leads to differential rotation in the interior in the time scales of order of $\sim\gamma^{-1}$. In this paper, we assumed that the star is uniformly rotating when computing tidal responses and that the time dependence of tidally driven mean flows is given by $e^{\gamma t}$ and that of the responses by $e^{\gamma t/2}$ to derive the governing equations for mean flows of second order. Probably, this is not necessarily a good approximation for the problem when we consider binary evolution in the time scales longer than $\sim\gamma^{-1}$, in which time scales equilibrium rotation laws would become substantially different from uniform rotation. It is thus highly desirable to follow time development of mean flows as a result of interactions between the mean flows and tidal responses in differentially rotating stars.

We should be cautious about the results suggested by Figure 12, where we have computed $v_{\rm max}^{(2)}/(R\sigma_{0})$ at $x=0.99$ as a function of $\bar{\omega}_{\rm tide}$ assuming that $\bar{\gamma}$ is a constant. For example, however, if we assume that $\bar{\gamma}$ is given by $\bar{\gamma}_{\rm tide}\propto{\cal T}_{2}$, $\bar{\gamma}$ will make sharp resonance peaks as a function of $\bar{\omega}_{\rm tide}$. The rapid increase in $|\bar{\gamma}|$ at peaks, on the other hand, will suppress the resonance peaks of $v_{\rm max}^{(2)}/(R\sigma_{0})$ when $v_{\rm max}^{(2)}\bar{\gamma}\sim C$ holds where $C$ is a constant that does not depend on $\bar{\omega}_{\rm tide}$. This suggests that the magnitudes of $v_{\rm max}^{(2)}$ will be only weakly dependent on $\bar{\omega}_{\rm tide}$ even near resonance although the flow patterns in resonance will be different from those in off-resonance. If this is the case for $\bar{\gamma}=\bar{\gamma}_{\rm tide}$, we can estimate the magnitudes of $v_{\rm max}^{(2)}/R\sigma_{0}$ at $x=0.99$ using the numerical results obtained for tidal resonance with low radial order $g$-modes. As suggested in the last paragraph of the previous section, for the parameter values of $q\sim 0.1$ and $a_{\rm orb}/R\sim 10$, for example, we have $\bar{\gamma}_{\rm tide}\sim 10^{-8}$, for which the magnitude of $v_{\rm max}^{(2)}/R\sigma_{0}$ will be of order of $\sim 10^{4}$. This value is too large to be accepted. Obviously we need more careful analyses concerning possible amplitudes of the mean flows driven by tidal responses.

In this paper, we assumed that the star is slowly rotating at $\bar{\Omega}=0.1$. For slow rotation, low radial order $g$-modes are not necessarily significantly affected by rotation, and there arise no significant differences in the mode properties between prograde and retrograde low radial order $g$-modes, although there appear on the retrograde side sequences of $r$-modes whose oscillation frequency in the co-rotating frame of the star is comparable to or less than $\Omega$. For rapidly rotating stars, as suggested by Figure 2, the tidal responses will have properties qualitatively different from those in slowly rotating stars, even in the frequency ranges of low radial order $g$-modes. The properties of tidal responses of a massive star also depend on the evolutional stages. As the star evolves from the ZAMS stage, the frequency spectra of low frequency $g$-modes will be denser and the amplitudes of $g$-modes tend to be confined into the deep interior. The development of a $\mu$-gradient zone outside the convective core will make the frequency spectra more complicated. Since low frequency $g$-modes can be trapped in the well-developed $\mu$-gradient zone, if the tidal forcing is in resonance with $g$-modes trapped in the $\mu$-zone, mean flows driven by the $g$-modes will have mixing effects on material there even if non-adiabatic effects are small in the deep interior.

Tidal responses and tidally driven mean flows of the star discussed in this paper have rather simple properties since no low frequency modes of the model are pulsationally unstable. Probably, this is not the case for slowly pulsating (SPB) stars, because numerous low frequency $g$-modes and $r$-modes of the stars are destabilized by the opacity bump mechanism. For these variable stars, there exists a strong excitation zone that surpasses damping contributions in the interior. The sign of the tidal torque ${\cal T}_{2}$ may change as a function of $\omega_{\rm tide}$. We expect that tidal mean flows driven in SPB stars will have different properties from those in massive main sequence stars.

Appendix A differential equations for tidal responses

Substituting the series expansions (25) to (28) into the perturbed basic equations (19), (20), (21), and (22), we obtain a finite set of linear ordinary differential equations for the expansion coefficients (see, e.g., Lee & Saio 1987). If we use vector notation for the set of differential equations, defining the dependent variables \mbox{\boldmathy}_{j}, $h$, $t$, and $\psi$ as

[TABLE]

we write the set of linear ordinary differential equations for tidally forced non-adiabatic oscillations of rotating stars as

[TABLE]

[TABLE]

[TABLE]

[TABLE]

where

[TABLE]

and

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

Note that the relations between the variables (\mbox{\boldmathh},{\rm i}\mbox{\boldmatht}) and (\mbox{\boldmathy}_{1},\mbox{\boldmathY}_{2}) are given

[TABLE]

[TABLE]

where

[TABLE]

and \mbox{\boldmath\sf O}^{T} is the transpose matrix of $\sf O$, $\sf I$ is the unit matrix. The non-zero elements of the matrices \mbox{\boldmath\sf\Lambda}_{0}, \mbox{\boldmath\sf\Lambda}_{1}, \mbox{\boldmath\sf L}_{0}, \mbox{\boldmath\sf L}_{1}, \mbox{\boldmath\sf M}_{0}, \mbox{\boldmath\sf M}_{1}, $\sf K$, \mbox{\boldmath\sf C}_{0} for even modes are

[TABLE]

[TABLE]

[TABLE]

and for odd modes

[TABLE]

[TABLE]

[TABLE]

where $l_{j}=2(j-1)+|m|$ and $l_{j}^{\prime}=l_{j}+1$ for even modes and $l_{j}=2j-1+|m|$ and $l_{j}^{\prime}=l_{j}-1$ for odd modes, and

[TABLE]

for $l\geq|m|$, and $J^{m}_{l}=0$ otherwise.

We note that the terms proportional to $\psi$ are inhomogeneous terms of the set of linear differential equations, and if we drop these inhomogeneous terms the set of linear ordinary differential equations reduce to those for free oscillations of stars (Lee & Saio 1987). The oscillation frequency $\omega$ should be regarded as the tidal forcing frequency $\omega_{\rm tide}$ for tidal responses.

To integrate the set of linear ordinary differential equations, we employ a Henyey type method of integration. For free oscillations of stars, for example, we formally write the set of linear differential equations as

[TABLE]

where $x=\ln r$,

[TABLE]

and $\sf C$ is the coefficient matrix. The differential equation (95) may reduce to a set of difference equations given by

[TABLE]

where $n$ is the mesh number of the background model, running from $n=1$ (the center) to $n=N$ (the surface of the model), and we usually assume $\alpha=1/2$. Equations (97) give recurrence equations

[TABLE]

where

[TABLE]

The inner and outer boundary conditions and the amplitude normalization may be written as

[TABLE]

where \mbox{\boldmath\sf B}_{\rm in} and \mbox{\boldmath\sf B}_{\rm out} are the coefficient matrices defining the boundary conditions. Using Newton-Raphson method, we look for $\omega$ such that the functions \mbox{\boldmathY}^{n} satisfy all the recurrence relations (98), the boundary conditions, and amplitude normalization (100). The background models we use in this paper have more than 2000 mesh points in the radial direction, which makes it possible for us to get accurate eigenmodes even when the modes have radial nodes of the eigenfunctions as many as $\sim 100$. Note that for tidally forced oscillations the vectors \mbox{\boldmathd}^{n} become nonzero vectors because of the inhomogeneous terms due to $\psi$ and that we omit the normalization $S_{l_{1}}^{N}=1$ to calculate forced oscillations.