© ESO, 2016

1. Introduction

Peters (1989) defines an Algol as a semidetached binary in which the less massive star fills its Roche lobe (RL) and is the cooler, fainter and larger, while the most massive star does not fill its RL and is still on the main sequence.

The evolution of a binary is conservative when mass and angular momentum of the system remain constant. Eggleton (2000) introduced the word “liberal” when this is not the case. However, it is difficult to calculate the amount of mass that leaves the system. Wellstein (2001) assumes that mass is lost from a massive system as soon as the gainer spins at critical velocity. Van Rensbergen et al. (2008) expand this criterion estimating the temperature of a high temperature accretion region (HTAR) so that the sum of rapid rotation and radiation pressure from the hot region overcomes gravitation and drives mass out of the binary. This criterion should be refined as soon as more observed values of sizes and temperatures of hot regions are available.

Orbital periods of Algols are known accurately. Somewhat different values of masses, radii and effective temperatures of the components are found in the literature (e.g. in the catalogs of Budding et al. 2004 and Brancewicz & Dworak 1980).

The most important problems in modeling an individual system are that a binary evolutionary calculation starts with assumed values of initial mass of the donor (subscript d) and the gainer (subscript g) and the initial period of the system. This introduces an intrinsic uncertainty since the present values of the masses are not accurately known. Very little to nothing is known about possible present or past mass loss and subsequent loss of angular momentum from the system.

Algol systems with an accretion disk are a special subset, as the accretion disk puts additional constraints on modeling the evolution of the system. The characteristics of six of them are elaborated further in Sect. 4. We will use them to constrain the parameters of the progenitors and their evolution. We first discuss the changes that have been made to the evolution code since the one used in previous liberal binary calculations (Van Rensbergen et al. 2008, 2011) and the joined catalog.

2. Modifications in the binary evolution code

2.1. Tides

Tidal interactions modulate the angular velocity of a binary member (ω) with the angular velocity of the system (ω_orb). Following Wellstein (2001) we have $\frac{\mathrm{1}}{\mathit{\omega }\mathrm{-}{\mathit{\omega }}_{\mathrm{orb}}}\frac{\mathrm{d}{\mathit{\omega }}_{\mathrm{orb}}}{\mathrm{d}\mathit{t}}\mathrm{=}\mathrm{-}\frac{\mathrm{1}}{{\mathit{\tau }}_{\mathrm{sync}}{\mathit{f}}_{\mathrm{sync}}}\mathrm{=}\mathrm{-}\frac{\mathrm{1}}{{\mathit{t}}_{\mathrm{sync}}}$(1)where t_sync is the synchronization time and τ_sync is a characteristic time associated with the theory of Darwin (1879) ${\mathit{\tau }}_{\mathrm{sync}}\mathrm{\left(}\mathrm{yr}\mathrm{\right)}\mathrm{=}{\mathit{q}}^{-2}{\left(\frac{\mathit{a}}{\mathit{R}}\right)}^{\mathrm{6}}\mathit{.}$(2)The star that is synchronized is in the denominator of the expression for the mass ratio q. Wellstein (2001) uses f_sync = 1 for weak tides and f_sync = 0.1 for strong tides in his scenario for liberal evolution of massive close binaries. Synchronization is strong on a star with a radius (R) comparable to the size of the semimajor axis of the binary (a).

If stars rotate asynchronously with the orbit, their spin changes by the tides, changing their angular momentum with an amount $\mathrm{\Delta }{\mathit{J}}_{\mathrm{spin}}\mathrm{=}\mathit{I}\mathrm{(}{\mathit{\omega }}_{\mathrm{orb}}\mathrm{-}\mathit{\omega }\mathrm{)}\left[\mathrm{1}\mathrm{-}{\mathrm{e}}^{\frac{\mathrm{-}\mathrm{\Delta }\mathit{t}}{{\mathit{t}}_{\mathrm{sync}}}}\right]\mathit{.}$(3)To avoid the artificial quantity f_sync we calculate t_sync with a physical model, discriminating between convective and radiative envelopes. Hurley et al. (2002) use an expression of Hut (1981) to calculate t_sync for circular orbits of stars with a convective envelope ${\mathit{t}}_{\mathrm{sync}}\mathrm{=}{\mathit{q}}^{-2}{\left(\frac{\mathit{a}}{\mathit{R}}\right)}^{\mathrm{6}}\frac{\mathrm{1}}{\mathrm{3}{\mathit{k}}_{\mathrm{2}}}{\mathit{t}}_{\mathrm{F}}\frac{\mathit{I}}{\mathit{M}{\mathit{R}}^{\mathrm{2}}}$(4)where k₂ is the apsidal motion constant; t_F is the viscous friction time; and M, R, and I are respectively the stellar mass, radius and moment of inertia. Following Zahn (1977) it is possible to write ${\mathit{t}}_{\mathrm{F}}\mathrm{=}\frac{\mathit{R}}{{\mathit{\nu }}_{\mathit{t}}}$, where ν_t is the coefficient of eddy viscosity. It is clear that ν_t = 0 for a star with an envelope in radiative equilibrium, so that for these stars t_F is infinite and tidal friction is not at work. If a convection region is a substantial fraction of a star and if convection transports most of the energy flux, the viscous friction time (Zahn 1977) is ${\mathit{t}}_{\mathrm{F}}\mathrm{\left(}\mathrm{yr}\mathrm{\right)}\mathrm{=}{\left(\frac{\mathit{M}{\mathit{R}}^{\mathrm{2}}}{\mathit{L}}\right)}^{\frac{\mathrm{1}}{\mathrm{3}}}\mathrm{=}\mathrm{0.4311}{\left(\frac{\frac{\mathit{M}}{{\mathit{M}}_{\mathrm{\odot }}}{\left[\frac{\mathit{R}}{{\mathit{R}}_{\mathrm{\odot }}}\right]}^{\mathrm{2}}}{\frac{\mathit{L}}{{\mathit{L}}_{\mathrm{\odot }}}}\right)}^{\frac{\mathrm{1}}{\mathrm{3}}}\mathrm{·}$(5)Relation (4) can thus be written as: ${\mathit{t}}_{\mathrm{sync}}\mathrm{\left(}\mathrm{yr}\mathrm{\right)}\mathrm{=}\mathrm{0.4311}{\mathit{q}}^{-2}{\left(\frac{\mathit{a}}{\mathit{R}}\right)}^{\mathrm{6}}\frac{\mathrm{1}}{\mathrm{3}{\mathit{k}}_{\mathrm{2}}}\frac{\mathit{I}}{\mathit{M}{\mathit{R}}^{\mathrm{2}}}{\left(\frac{\frac{\mathit{M}}{{\mathit{M}}_{\mathrm{\odot }}}{\left[\frac{\mathit{R}}{{\mathit{R}}_{\mathrm{\odot }}}\right]}^{\mathrm{2}}}{\frac{\mathit{L}}{{\mathit{L}}_{\mathrm{\odot }}}}\right)}^{\frac{\mathrm{1}}{\mathrm{3}}}\mathrm{·}$(6)During stellar evolution k₂ can be calculated from a simplified formula given by Kopal (1959), also used by Odell (1974) and Kopal (1978), using an integral from the center to the edge of the gainer. ${\mathit{k}}_{\mathrm{2}}\mathrm{=}\frac{\mathrm{16}\mathit{\pi }}{\mathrm{5}\mathit{M}{\mathit{R}}^{\mathrm{2}}}{\mathrm{\int }}_{\mathrm{0}}^{{\mathit{R}}_{\mathrm{g}}}\mathit{\rho }\mathrm{\left(}\mathit{r}\mathrm{\right)}{\mathit{r}}^{\mathrm{2}}\mathrm{d}\mathit{r}\mathit{.}$(7)Hilditch (2001) gives an expression of Zahn (1975) to calculate the tidal action on a star with a radiative envelope and a convective core. Owing to the convective core the tidal torque mechanism is at work resulting from the oscillation of the convective core. The corresponding synchronization time is $\begin{array}{ccc}{\mathit{t}}_{\mathrm{sync}}\mathrm{\left(}\mathrm{yr}\mathrm{\right)}& \mathrm{=}& \mathrm{3.1816}\mathrm{\times }{\mathrm{10}}^{-6}{\left(\frac{\frac{\mathit{R}}{{\mathit{R}}_{\mathrm{\odot }}}}{\frac{\mathit{M}}{{\mathit{M}}_{\mathrm{\odot }}}}\right)}^{\frac{\mathrm{1}}{\mathrm{2}}}\frac{\mathit{I}}{\mathit{M}{\mathit{R}}^{\mathrm{2}}}\\ & & \end{array}$(8)The strength of the tide depends critically upon the size of the convective core. The tidal-torque constant E₂ is zero for stars with mass below 1.25 M_⊙ so that this mechanism has no effect on a star with a convective atmosphere.

However, the radius of the convective core is larger than given by the Schwarzschild (1958) criterion. Using the mixing length ℓ of the path after which a convective cell dissolves and the pressure scale height H_P, the overshooting parameter α is defined as α = ℓ/H_P.

Our choice of values of E₂ calculated by Claret (2004) with overshooting of the convective core characterized by α = 0.2 is justified for the sample of binaries in this paper. We used the spectroscopic Hertzsprung-Russell diagram (HRD) proposed by Langer & Kudritzki (2014). The location in this diagram of the “red points” (coolest points during main sequence evolution) of massive stars was taken from Castro et al. (2014) and MacDonald et al. (2012) for lower mass stars. We calculated the main sequence evolution of single stars between 2 and 30 M_⊙ using different values of the parameter α. We found that the choice α = 0.2 is justified for stars with masses below 12.5 M_⊙, growing gradually up to α = 0.4 at 30 M_⊙.

Expressions (6) and (8) refine and quantify the original theory of Darwin (1879). Belczynski et al. (2008) apply radiative damping to stars with radiative envelopes: Main sequence stars with mass above M_MS,conv = 1.15 M_⊙, CHeB stars above 7 M_⊙, massive evolved He stars, and He MS stars. For all other stars convective damping is applied. We assume that during Roche lobe overflow (RLOF) the atmospheres of both stars are sufficiently disturbed so they can be considered to be in convective equilibrium.

Tidal action produced by meridional circulation (Tassoul 2000) can be added to the Darwin theory in order to realize synchronization of binaries over a wider range than performed by the Darwin interaction alone. The synchronization time in the Tassoul theory is given by $\begin{array}{ccc}{\mathit{t}}_{\mathrm{sync}}\mathrm{\left(}\mathrm{yr}\mathrm{\right)}& \mathrm{=}& \mathrm{5.396}\mathrm{\times }{\mathrm{10}}^{\mathrm{2}\mathrm{+}\mathit{\sigma }\mathrm{-}\frac{\mathit{N}}{\mathrm{4}}}\frac{\mathrm{1}\mathrm{+}\mathit{q}}{\mathit{q}}{\left(\frac{{\mathit{L}}_{\mathrm{\odot }}}{\mathit{L}}\right)}^{\frac{\mathrm{1}}{\mathrm{4}}}{\left(\frac{\mathit{M}}{{\mathit{M}}_{\mathrm{\odot }}}\right)}^{\frac{\mathrm{3}}{\mathrm{4}}}\\ & & \end{array}$(9)In the same way as t_F is the main uncertainty of the Hut theory in Eq. (4), the values of N and σ are the main uncertainties of the Tassoul theory in Eq. (9). The number N is given as a function of the coefficients of turbulent and radiative viscosity: ${\mathrm{10}}^{\mathit{N}}\mathrm{=}\frac{{\mathit{\nu }}_{\mathrm{rad}}\mathrm{+}{\mathit{\nu }}_{\mathrm{turb}}}{{\mathit{\nu }}_{\mathrm{rad}}}$. In the absence of turbulence we have for a radiative envelope N = 0. We assume with Tassoul & Tassoul (1992) that N takes the rather large value N = 10 for a convective envelope, because turbulent viscosity is always much larger than radiative viscosity if convection is at work. The always uncertain value of σ will be calculated within these assumptions so that relation (9) reproduces observations best.

Using Kepler’s 3rd law t_sync increases in relation (9) with increasing separations as ${{}^{\mathrm{\right(}}\frac{\mathit{a}}{\mathit{R}}^{\mathrm{\left)}}}^{\mathrm{4.125}}$ which is much slower than ${{}^{\mathrm{\right(}}\frac{\mathit{a}}{\mathit{R}}^{\mathrm{\left)}}}^{\mathrm{6}}$ in relation (6) and ${{}^{\mathrm{\right(}}\frac{\mathit{a}}{\mathit{R}}^{\mathrm{\left)}}}^{\mathrm{8.5}}$ in relation (8).

Claret et al. (1995) propose σ = 1.6, leading to values of the exponent respectively 3.6 (radiative) and 1.1 (convective) in relation (9). We propose calculating σ so as to keep the components of binaries synchronized up to orbital periods of 10 days as observed by Matthews & Mathieu (1992) for main sequence A stars. We calculated an extended set of binary main sequence evolutions for stars between 1 and 30 M_⊙ with lower mass companions. Initial periods were between 5 and 15 days, so that no RLOF occurred during the main sequence life of the most massive component. The binaries exerted tidal interaction on one another, starting with the Darwin tide (relation (6) for stars with a convective envelope and relation (8) for stars with a radiative envelope). The Darwin mechanism did not keep the observed fraction of binaries synchronized. We added the Tassoul tide (9) using different values of σ in order to achieve synchronization within reasonable limits.

Figure 1 shows the evolution of the factor SYNC = $\frac{{\mathit{\omega }}_{\mathrm{spin}}}{{\mathit{\omega }}_{\mathrm{orb}}}$ for the primary in a binary pair 5 M_⊙ + 2 M_⊙. The system starts synchronized and the evolution is followed from zero age main sequence (ZAMS) to the end of hydrogen core burning for the primary. During this time the less massive component remains synchronized, even if no tides are at work. Without tide, the 5 M_⊙ star is soon asynchronized. Accurate inspection of Fig. 1 shows that the Darwin tide contributes almost nothing to the synchronization of the 5 M_⊙ primary. The contribution of the weak Tassoul tide (σ = 4.5) is not of very much help. The Darwin tide enforced by the strong Tassoul tide (σ = 3.5) restores synchronization most of the time for this system with an orbital period of 10 days. Therefore, the value $\mathrm{2}\mathrm{+}\mathit{\sigma }\mathrm{-}\frac{\mathit{N}}{\mathrm{4}}\mathrm{=}\mathrm{5.5}$ will be used in relation (9) for the tidal tuning of a star with a radiative atmosphere.

Fig. 1 Evolution of SYNC = $\frac{{\mathit{\omega }}_{\mathrm{spin}}}{{\mathit{\omega }}_{\mathrm{orb}}}$, during the main sequence life of the 5 M⊙ primary in a 5 M⊙ + 2 M⊙ binary with a 10-day orbital period. From bottom to top: without tide, Darwin tide contributing almost nothing to the synchronization, weak Tassoul tide added, strong Tassoul tide added. The 5 M⊙ gainer needs the strong Tassoul tide to remain sufficiently synchronized.

In order to find a value for σ for a star with a convective atmosphere we followed the evolution of a binary pair 1.15 M_⊙ + 1 M_⊙. We find that without tide the 1.15 M_⊙ star is soon asynchronized. Convective Darwin tides restore synchronization for a 10 d orbital period without help of the Tassoul tide. However, for a 15 d orbital period, the Darwin tide restores synchronization less well. The Tassoul tide in its weak form (σ = 5) restores synchronization completely. Therefore, we use the value $\mathrm{2}\mathrm{+}\mathit{\sigma }\mathrm{-}\frac{\mathit{N}}{\mathrm{4}}\mathrm{=}\mathrm{4.5}$ in relation (9) for the tidal tuning of stars with a convective atmosphere.

The Tassoul and Darwin mechanisms act simultaneously, so that the resulting synchronization times are given by $\frac{\mathrm{1}}{{\mathit{t}}_{\mathrm{sync}}}\mathrm{=}\frac{\mathrm{1}}{{\mathit{t}}_{\mathrm{sync}\mathit{,}\mathrm{Darwin}}}\mathrm{+}\frac{\mathrm{1}}{{\mathit{t}}_{\mathrm{sync}\mathit{,}\mathrm{Tassoul}}}\mathrm{·}$(10)

2.2. Conservation of angular momentum

The total angular momentum J_Σ of a binary is the sum of the orbital angular momentum J_orb and the spin angular momenta of gainer J_g and donor J_d. Tides continuously exchange amounts of ΔJ_d , ΔJ_g ,ΔJ_orb. In the conservative case $\mathrm{\Delta }{\mathit{J}}_{\mathrm{\Sigma }}\mathrm{=}\mathrm{\Delta }{\mathit{J}}_{\mathrm{d}}\mathrm{+}\mathrm{\Delta }{\mathit{J}}_{\mathrm{g}}\mathrm{+}\mathrm{\Delta }{\mathit{J}}_{\mathrm{orb}}\mathrm{=}\mathrm{0.}$(11)However, the system loses angular momentum through stellar wind (ΔJ_wind< 0), and during its liberal era also through mass loss out of the system (ΔJ_out< 0). Conservation of angular momentum during evolution is then imposed by $\mathrm{\Delta }{\mathit{J}}_{\mathrm{d}}\mathrm{+}\mathrm{\Delta }{\mathit{J}}_{\mathrm{g}}\mathrm{+}\mathrm{\Delta }{\mathit{J}}_{\mathrm{orb}}\mathrm{-}\mathrm{\Delta }{\mathit{J}}_{\mathrm{wind}}\mathrm{-}\mathrm{\Delta }{\mathit{J}}_{\mathrm{out}}\mathrm{=}\mathrm{0}$(12)In this paper, we calculated ΔJ_wind using Vink et al. (2001) for stars hotter than 12500 K and De Jager et al. (1988) for cooler stars. In the case of liberal evolution the quantity ΔJ_out was calculated as in previous papers (Van Rensbergen et al. 2008, 2011; see also Siess et al. 2013) assuming that mass lost from the system takes only the specific orbital angular momentum of the gainer.

3. Observations

The observed data for which we calculated the most plausible progenitors are given in Tables 1 (system parameters) and 2 (disk parameters). Table 1 gives the observed system parameters to be met by the evolution of our progenitors. Masses and radii of β Lyr are from Zhao et al. (2008), the effective temperatures and mass transfer rates from Harmanec (2002). Masses, radii and effective temperatures for AU Mon are from Desmet et al. (2010) and Atwood-Stone et al. (2012). We determined a very uncertain mass transfer rate applying the (O−C) method on the data found in the catalog of Kreiner et al. (2001). Masses, radii and effective temperatures for V356 Sgr are from Dominis et al. (2005). The mass transfer rate is from Ziolkowski (1985). Masses, radii and effective temperatures for TT Hya are from Miller et al. (2007). The mass transfer rate was determined from the data found in the catalog of Kreiner et al. (2001). Masses, radii and effective temperatures for RY Per are from Peters & Polidan (2004). No trustworthy determination of the mass transfer rate was found in this case. Masses and radii for SW Cyg are from Richards & Albright (1999), the effective temperatures from Budding et al. (2004) whereas the mass transfer rate is from Qian et al. (2002).

Table 2 compares observed luminosities, temperatures and sizes of the six disks with the same quantities obtained in Sect. 4.2 from the best models. Disk luminosities for SW Cyg and TT Hya are taken from Albright & Richards (1996). The outer edge temperature of AU Mon is from Djurasevic et al. (2010). Other characteristic disk temperatures are from Miller et al. (2007) for TT Hya, from Harmanec (2002) and Ak et al. (2007) for β Lyr and from Wilson & Caldwell (1978) for V356 Sgr. The sizes are from Atwood-Stone et al. (2012) for AU Mon, Albright & Richards (1996) for SW Cyg, Peters (2001) for TT Hya, Sudar et al. (2011) for RY Per, Mennickent & Djurasevic (2013) for β Lyr and from Wilson & Caldwell (1978) for V356 Sgr. The size of the disk in Table 2 is the disk radius divided by the Roche radius of the gainer. Since accretion disks in Algols are seen almost edge on, the model temperatures that were calculated at the edge show a fair agreement with the observations in Table 2.

4. Modeling six Algols with accretion disks

4.1. Modeling the stars

Using our binary evolutionary code with the modifications discussed in Sect. 2 we determine in this subsection the progenitors of six binaries with an accretion disk around the gainer. Four of them turn out to evolve conservatively (AU Mon, SW Cyg, TT Hya and RY Per), whereas β Lyr and V356 Sgr have liberal eras during their evolution. In the conservative case the initial orbital period of a system follows directly from the present situation (subscript 1) for any combination of the initial masses (subscript 0). $\frac{{\mathit{P}}_{\mathrm{1}}}{{\mathit{P}}_{\mathrm{0}}}\mathrm{=}{\left(\frac{{\mathit{M}}_{\mathrm{d}\mathit{,}\mathrm{0}}{\mathit{M}}_{\mathrm{g}\mathit{,}\mathrm{0}}}{{\mathit{M}}_{\mathrm{d}\mathit{,}\mathrm{1}}{\mathit{M}}_{\mathrm{g}\mathit{,}\mathrm{1}}}\right)}^{\mathrm{3}}\mathrm{·}$(13)This relation disregards small deviations that are caused by the stellar wind of both components and the action of the tides. Mennickent (2014) proposes a (4 M_⊙ + 3.6 M_⊙, P = 3 d) progenitor for AU Mon with a 1.533 M_⊙ donor and a 6.067 M_⊙ gainer at present. Since we model AU Mon towards a 1.2 M_⊙ donor and a 7 M_⊙ gainer our progenitors with initial mass ratio (4/3.6) need to have an initial orbital period a little below 1.4 d, determined from expression (13). With this small initial orbital period the stars in the binary merge soon after birth.

All the gainers in Table 1 are main sequence stars whereas donors are evolved beyond the main sequence.

Initial orbital periods of conservative cases were calculated with expression (13).

Loss of mass in the liberal cases is triggered by the combined action of rapid rotation and radiation pressure from the HTAR region where the transformation of matter energy into radiation is crucial. The radiative efficiency is defined through the relation: L_add = ηṀc². For comparison we give some numbers: η ≈ 0.007 for nuclear fusion, η ≈ 10^-4 for an accreting white dwarf, η ≈ 10^-2 for an accreting neutron star. The value of η for a hot spot on a main sequence gainer in the case of direct impact was calculated by Van Rensbergen et al. (2011) and found to be significantly smaller than in the case of a mass accreting white dwarf. In the case of a gainer surrounded by a disk, the hot spot could be located at the edge of the disk. Moreover, Bisikalo (2007) observed and reproduced hot lines using 2D and 3D gas dynamical simulations for disk systems. Using the same criterion as used by Van Rensbergen et al. (2011) we found that four out of six systems treated in this paper evolved conservatively. The two liberal systems (V356 Sgr and β Lyr) were observed to be ejecting mass into space.

Initial periods for liberal cases were calculated with Eq. (12) assuming a loss of specific angular momentum ΔJ_out equal to the orbital angular momentum of the gainer. The mass loss out of the system was found to be vertical on the orbital plane and has been interpreted as bipolar jets for V356 Sgr by Peters & Polidan (2004) and for β Lyr by Ak et al. (2007) and Harmanec (2002). Mass loss through the second Lagrangian point is unlikely because this point is far away from the outer edge of the accretion disk, and the large amount of angular momentum lost narrows the orbit so that the formation of an accretion disk around the gainer is inhibited.

Liberal cases were calculated with various amounts of mass lost by the system. However, assuming a mass of [1.5−2] M_⊙ lost by β Lyr and V356 Sgr gives the best fit with the observed data.

The progenitor systems are characterized by (P₀,M_d,0,M_g,0). Evolutionary sequences were calculated for various combinations of (M_d,0,M_g,0) and evaluated at P = P_obs. The results are shown in Table 3.

All the binaries in Table 1 have a mass ratio q< 0.27. For those binaries that started RLOF during core hydrogen burning of the donor this means that they are at the end of their RLOF evolution. TT Hya and RY Per are not in that case.

4.2. Modeling the disk

Our code also reproduces characteristics of the accretion disk around the gainer. The criterion of Lubow & Shu (1975) uses the concept of fractional radius ϖ which is the radius divided by the semimajor axis. The fractional radius is determined at the edge of the gainer ϖ_g and at specific distances ϖ_d¹ and ϖ_min so that a semidetached binary changes from direct hit into accretion through a transient disk as soon as ϖ_g = ϖ_d and later on through a permanent accretion disk when ϖ_g = ϖ_min. At ϖ_d the disk radius starts from R_disk = R_g until the radius of a permanent disk R_PD is achieved at ϖ_min. For the radius of a permanent disk we assumed, as do Carroll & Ostlie (2007), that R_PD = 2 r_circ, a value that is a somewhat larger than 2 ϖ_d.

The radius of the transient disk is calculated with the assumption that the disk radius grows faster than linearly from ϖ_d and R_g towards ϖ_min and R_PD. A slower growth yields sizes of transient disks that are smaller than observed.

The Lubow & Shu Fig. 2 shows the fractional radii of the six Algol gainers in this paper. All the observed fractional radii of the gainers are located below ϖ_d as needed to develop an accretion disk. TT Hya is located below ϖ_min. The calculated disk radii are listed in Table 2. The gainers of four systems are surrounded by a transient accretion disk. TT Hya has a permanent disk and β Lyr is at the transition between the transient and permanent disk. The observed size of the accretion disk in SW Cyg is much larger than the calculated one. Since the present position of SW Cyg in Fig. 2 suggests that SW Cyg is far from the permanent disk phase, we conclude that the observed size overestimates reality.

Fig. 2 Fractional radii (radius/semimajor axis) of six gainers and their accretion disks in the Lubow & Shu diagram. The observed radii of gainers are dots below ϖd. The calculated disk radii are squares above ϖd. The two observed disk radii are triangles. The fractional radii of permanent disks and the position of the gainer’s Roche radius are also indicated.

The temperature in the disk is both caused by accretion and by the luminosity of the gainer. The accretion luminosity of a disk can be approximated with $\frac{{\mathit{L}}_{\mathrm{acc}\mathit{,}\mathrm{disk}}}{{\mathit{L}}_{\mathrm{\odot }}}\mathrm{=}\frac{\mathit{G}{\mathit{M}}_{\mathrm{g}}\mathit{Ṁ}\mathrm{d}}{\mathrm{2}{\mathit{R}}_{\mathrm{g}}}\mathrm{=}\mathrm{1.57066}\mathrm{\times }{\mathrm{10}}^{\mathrm{7}}\frac{\left(\frac{{\mathit{M}}_{\mathrm{g}}}{{\mathit{M}}_{\mathrm{\odot }}}\right)\mathit{Ṁ}\mathrm{d}\left(\mathrm{in}\frac{{\mathit{M}}_{\mathrm{\odot }}}{\mathrm{yr}}\right)}{\left(\frac{{\mathit{R}}_{\mathrm{g}}}{{\mathit{R}}_{\mathrm{\odot }}}\right)}\mathrm{·}$(14)This expression can be calculated even if the gainer has no disk. Hence, it is better to calculate the accretion luminosity considering Bath & Pringle (1981), Horne (1985) and Rutten et al. (1992) who give the distribution of the accretion temperature through a disk around a white dwarf or a neutron star, using a characteristic disk temperature T_disk. With masses and radii in solar units and mass transfer rates in M_⊙/yr, this characteristic disk temperature created by accretion is defined by ${\mathit{T}}_{\mathrm{disk}}\mathrm{=}\mathrm{478074}{\left(\frac{{\mathit{M}}_{\mathrm{g}}\mathit{Ṁ}\mathrm{d}}{{\mathit{R}}_{\mathrm{g}}^{\mathrm{3}}}\right)}^{\mathrm{0.25}}\mathrm{·}$(15)The distribution of the temperature through the disk is given by ${\mathit{T}\mathrm{\left(}\mathit{r}\mathrm{\right)}}^{\mathrm{4}}\mathrm{=}{\mathit{T}}_{\mathrm{disk}}^{\mathrm{4}}{\left(\frac{{\mathit{R}}_{\mathrm{g}}}{\mathit{r}}\right)}^{\mathrm{3}}\left(\mathrm{1}\mathrm{-}\sqrt{\frac{{\mathit{R}}_{\mathrm{g}}}{\mathit{r}}}\right)\mathrm{·}$(16)Integrating over both sides of the disk, the amount of accretion luminosity in a ring with width dr at a distance r from the gainer is $\mathrm{d}\mathit{L}\mathrm{=}\mathrm{4}\mathit{\pi }{\mathit{\sigma }}_{\mathrm{R}}\mathit{r}\mathrm{d}\mathit{r}\mathit{T}\mathrm{(}\mathit{r}{\mathrm{)}}^{\mathrm{4}}\mathrm{erg}{\mathrm{s}}^{-1}\mathrm{·}$(17)Integrating this quantity from the gainer’s edge to the outer edge of the disk, one finds $\frac{{\mathit{L}}_{\mathrm{acc}}}{{\mathit{L}}_{\mathrm{\odot }}}\mathrm{=}\mathrm{9.0203}\mathrm{\times }{\mathrm{10}}^{-16}{\mathit{T}}_{\mathrm{disk}}^{\mathrm{4}}{\left(\frac{{\mathit{R}}_{\mathrm{g}}}{{\mathit{R}}_{\mathrm{\odot }}}\right)}^{\mathrm{3}}\left(\frac{\mathrm{1}}{\mathrm{3}\frac{{\mathit{R}}_{\mathrm{g}}}{{\mathit{R}}_{\mathrm{\odot }}}}\mathrm{-}\frac{\mathrm{1}\mathrm{-}\frac{\mathrm{2}}{\mathrm{3}}\sqrt{\frac{{\mathit{R}}_{\mathrm{g}}}{{\mathit{R}}_{\mathrm{disk}}}}}{\frac{{\mathit{R}}_{\mathrm{disk}}}{{\mathit{R}}_{\mathrm{\odot }}}}\right)\mathrm{·}$(18)Restricting the luminosity of the disk to only accretion yields disk temperatures that are far below the temperatures listed in Table 2.

In Algols hot disks surround hot gainers, whereas cooler gainers yield cooler disks.

Moreover, the temperature at the inner edge of the disk needs to equal the effective temperature of the gainer, a requirement that is not fulfilled by relation (16). The temperature of the disk is thus also influenced by the outgoing radiation of the gainer. If the local luminosity in a disk is weakened only because of its distance from the source (i.e. there is no absorption), the temperature distribution is ${\mathit{T}}_{\mathrm{rad}\mathit{,}\mathrm{0}}\mathrm{=}{\mathit{T}}_{\mathrm{eff}}{\left(\frac{{\mathit{R}}_{\mathrm{g}}}{\mathit{r}}\right)}^{\mathrm{0.5}}\mathrm{·}$(19)Using this zero-order temperature distribution and integrating from the edge of the gainer to the outer edge of the disk, the radiation luminosity obtained is $\frac{{\mathit{L}}_{\mathrm{rad}}}{{\mathit{L}}_{\mathrm{\odot }}}\mathrm{=}\mathrm{9.0203}\mathrm{\times }{\mathrm{10}}^{-16}{\mathit{T}}_{\mathrm{eff}\mathit{,}\mathrm{g}}^{\mathrm{4}}{\left(\frac{{\mathit{R}}_{\mathrm{g}}}{{\mathit{R}}_{\mathrm{\odot }}}\right)}^{\mathrm{2}}\ln \frac{{\mathit{R}}_{\mathrm{disk}}}{{\mathit{R}}_{\mathrm{g}}}\mathrm{·}$(20)This expression overestimates the contribution of the radiation since disk matter will prevent outgoing radiation from traveling freely. Moreover, the temperatures measured at the outer edge of the disk are far below those calculated with this assumption. Unfortunately, an evaluation of the temperature at the outer edge has only been reported in the case of AU Mon (Djurasevic et al. 2010). Using their numbers the following is obtained $\frac{{\mathit{T}}_{\mathrm{rad}\mathit{,}\mathrm{edge}\mathit{,}\mathrm{real}}}{{\mathit{T}}_{\mathrm{rad}\mathit{,}\mathrm{edge}\mathit{,}\mathrm{0}}}\mathrm{=}\mathrm{0.6}\mathit{.}$(21)Since T_rad,edge,0 can easily be calculated with relation (19) we have adopted relation (21) to evaluate T_{rad,edge,real} for the five other systems. However, Table 2 shows only a moderate agreement between observed and calculated luminosities and temperatures of the disks. In order to improve the model, more and better determined values of these quantities are needed.

For our final model, we approximate the temperature distribution in the disk adding an exponentially decreasing function starting at the gainer’s surface (h = 0) where h is the distance above the gainer’s surface ${\mathit{T}}_{\mathrm{rad}}\mathrm{\left(}\mathit{h}\mathrm{\right)}\mathrm{=}{\mathit{T}}_{\mathrm{eff}\mathit{,}\mathrm{g}}{\left(\frac{{\mathit{R}}_{\mathrm{g}}}{\mathit{r}}\right)}^{\mathrm{0.5}}{\mathrm{e}}^{\frac{\mathrm{-}\mathit{h}}{\mathit{H}}}\mathrm{·}$(22)This temperature scale height H is determined so that expression (22) meets the effective temperature of the gainer at h = 0 and T_{rad,edge,real} at the outer edge of the disk $\mathit{H}\mathrm{=}\frac{\mathrm{-}\mathrm{(}{\mathit{R}}_{\mathrm{disk}}\mathrm{-}{\mathit{R}}_{\mathrm{g}}\mathrm{)}}{\ln \mathrm{0.6}}\mathrm{·}$(23)In this case, the integrated luminosity of the disk is $\frac{{\mathit{L}}_{\mathrm{rad}}}{{\mathit{L}}_{\mathrm{\odot }}}\mathrm{=}\mathrm{9.0203}\mathrm{\times }{\mathrm{10}}^{-16}{\mathit{T}}_{\mathrm{eff}\mathit{,}\mathrm{g}}^{\mathrm{4}}{\left(\frac{{\mathit{R}}_{\mathrm{g}}}{{\mathit{R}}_{\mathrm{\odot }}}\right)}^{\mathrm{2}}{\mathrm{e}}^{\frac{\mathrm{4}\mathit{R}\mathrm{g}}{\mathit{H}}}{\mathrm{\int }}_{{\mathit{R}}_{\mathrm{g}}}^{{\mathit{R}}_{\mathrm{disk}}}\frac{{\mathrm{e}}^{\frac{\mathrm{-}\mathrm{4}\mathit{r}}{\mathit{H}}}}{\mathit{r}}\mathrm{d}\mathit{r}\mathrm{·}$(24)The integral in relation (24) is the difference of the values of the exponential integral Ei($\frac{-4\mathit{r}}{\mathit{H}}$) at (r = R_disk) and (r = R_g). The radiative luminosity is added to the accretion luminosity (18) in order to obtain the real luminosity. The effective temperature of the disk, as given in Table 2, is calculated with ${\mathit{T}}_{\mathrm{eff}\mathit{,}\mathrm{disk}}\mathrm{=}\mathrm{5770}{\left(\frac{{\mathit{L}}_{\mathrm{disk}}}{{\mathit{L}}_{\mathrm{\odot }}}\right)}^{\frac{\mathrm{1}}{\mathrm{4}}}{\left(\frac{{\mathit{R}}_{\mathrm{\odot }}}{{\mathit{R}}_{\mathrm{disk}}}\right)}^{\frac{\mathrm{1}}{\mathrm{2}}}\mathrm{4}\sqrt{\mathrm{2}}\mathrm{·}$(25)

5. Evolution of the individual systems

Table 3 shows the characteristics of the initial binaries evolving best into the present characteristics given in Table 1. The evolutionary sequences possibly fitting the present observations were evaluated at the moment that the present orbital period was attained. The best set of initial parameters is always surrounded by a few different initial systems that evolve into current comparable results. However, progenitors with somewhat smaller initial mass of the donor do not develop the observed accretion disk, whereas future donor stars with a somewhat larger initial mass merge with the gainer before an accretion disk can be formed. In general we find that the masses of the initial system are defined within 0.2 M_⊙. We find that the initial periods are in the range [1.6−2.7] d, whereas the present orbital periods are in the range [4.5−13] d. The initial mass ratios are in the range q = M_d/M_g∈ [1.48−2.75] whereas the present mass ratios are in the range [0.17−0.27]. The present state of the binaries is thus far behind the orbital period minimum when both components of the binary have equal masses.

Table 4 illustrates the disk behavior and the life of the binary as an Algol. Table 4 mentions the duration of RLOF, and the percent of time that the disk appeared during RLOF. The percent of the total time of disk appearance that is shown by the present disk is also included. The binary system is always an Algol when the gainer is surrounded by an accretion disk. Four systems (AU Mon, β Lyr, V356 Sgr and SW Cyg) underwent previously RLOF during core hydrogen burning of the donor. Their disks live for a short fraction of the RLOF era. The other two systems (TY Hya and RY Per) start RLOF when hydrogen is exhausted in the core of the donor star. Their disks live for a large fraction of the following case B of RLOF. For all six cases, the present binary has a gainer still on the main sequence and a donor with a hydrogen exhausted core. The binary is an Algol during the largest part of the RLOF era. Inspection of Table 4 also shows that the four systems that lived their first RLOF during hydrogen core burning of the donor are currently (during their second RLOF) at the end of their evolution (β Lyr, AU Mon, V 356 Sgr and SW Cyg). The two systems that undergo their first and only RLOF during hydrogen shell burning of the donor are in the middle of their RLOF life (TT Hya and RY Per).

For each of the systems we calculated the evolution of the luminosity of the disk as a function of time. The accretion luminosity is calculated with relations (14) and (18). The evolution with time of the accretion luminosity reflects the evolution of the mass loss rate.

The radiation luminosity that contributes most to the luminosity of the disk is triggered by the radiation from the gainer that is absorbed in the disk and reradiated. Two approaches can be considered: the optically thin approach given by relation (20) and the optically thick approach given by relation (24). The resulting effect of the two approaches, added to the accretion luminosity (18), are shown in Figs. 3 and 4.

They show the calculated evolution of the luminosity of the disk divided by the luminosity of the gainer over the relative existence time of the disk, normalized to unity, respectively for the short lived disk and liberally evolving medium mass binary β Lyr and the longer lived disk and conservatively evolving low mass binary TT Hya. The vertical bar indicates the present disk with its present luminosity. It can be observed that in disks of Algol-type systems the accretion luminosity is of secondary importance as the luminosity originating from the gainer’s radiation is two to three orders of magnitude larger. The present luminosity of the disk of β Lyr (at 27% of the disk life time) is Log (L_disk/L_⊙) = 4.21, about 71% the luminosity of the gainer. When 48% of the disk lifetime has been consumed the luminosity of the disk with a radius of 35 R_⊙ will exceed the luminosity of the gainer. The present luminosity of the disk of TT Hya (at 21% of the disk life time) is L_disk/L_⊙ = 56, about 90% the luminosity of the gainer. When 30% of the disk lifetime has been consumed the luminosity of the large disk with a radius of 11.9 R_⊙ will exceed the luminosity of the gainer.

Fig. 3 Evolution of the luminosity of the short-lived disk of β Lyr divided by the gainer’s luminosity as a function of relative disk time. From bottom to top: accretion luminosities (14), (18), radiative luminosities (24) and (20) added.

Fig. 4 Evolution of the luminosity of the long-lived disk of TT Hya divided by the gainer’s luminosity as a function of relative disk time. From bottom to top: accretion luminosities (14), (18), radiative luminosities (24) and (20) added.

Fig. 5 Evolution of the gainer in β Lyr through the HR-diagram. The short-lived disk is the thick part on the evolutionary path. The best fit model is the dot on the thick path whereas the observed position of the gainer is the dot to the lower left.

Fig. 6 Evolution of the gainer of TT Hya. The dot on the thick path is the best model, whereas the observed position of the gainer is to the right.

Figures 5 and 6 show the evolution of the gainer in the HRD for the β Lyr and TT Hya. The disk phase is also indicated. β Lyr has lived a Case A era of direct impact on the gainer. TT Hya started RLOF after exhaustion of hydrogen in the core of the donor. Both donors evolved through the main sequence and are now transferring mass during their hydrogen shell burning. As the orbit of β Lyr is now rapidly widening a disk appears 2.63 × 10⁷ yr after ZAMS. This disk remains for 1.26 × 10⁵ yr, slightly more than 1% of the RLOF time. The system lives 9.551 × 10⁶ yr as an Algol. Previously published progenitors of β Lyr are from Ziolkowski (1985) (10 M_⊙ + 3.7 M_⊙, P = 3.44 d), De Greve & Linnell (1994) (9 M_⊙ + 7.65 M_⊙, P = 4 d) and Mennickent & Djurasevic (2013) (12 M_⊙ + 7.2 M_⊙, P = 2.5 d). The present initial masses approximate those of De Greve & Linnell (1994) calculated with external mass loss (ΔM = 0.45 M_⊙), but with a shorter initial period due to different values of external mass loss (ΔM = 1.5 M_⊙) and tidal interaction effects. The present model agrees well with the progenitor proposed by Mennickent & Djurasevic (2013), calculated with a somewhat larger amount of mass lost by the system.

In the case of TT Hya a disk appears 7.49 × 10⁸ yr after ZAMS. This disk remains for 1.60 × 10⁷ yr, approximately 87% of the RLOF time. The system lives 1.76 × 10⁷ yr as an Algol.

6. Conclusions

We extended our binary evolutionary code (Van Rensbergen et al. 2008, 2011) allowing for the action of viscous friction, tidal torque and meridional circulation. Matter leaving the system takes away the angular momentum of the gainer into space. Matter streams into interstellar space as a consequence of rapid rotation and the creation of high temperature accretion regions. The gauging of this event is still a matter of debate, depending on the well-known accretion luminosity and the poorly known size of the accretion region. Applying strictly the conservation of angular momentum we calculated the evolution of six progenitors towards six well-know Algols with accretion disks. We followed the evolution of a progenitor into the presently defined physical data of the binary components and the accretion disk around the gainer. The temperature of the disk is determined by accretion and by the radiation from the gainer’s surface. Fitting this rather large number of observed data allowed us to determine the masses and initial periods of the progenitors precisely. It has to be stated that different determinations of present mass, radii and effective temperatures of Algols would lead to different progenitors as found in this paper. Finally, we can summarize our findings as follows: We obtained a more accurate determination of initial parameters of observed Algol systems with a disk, by adding the disk parameters in the comparison with the models. Meridional circulation has to be included in the calculation of tidal interaction, at least for stars with radiative atmospheres. The luminosity of accretion disks is determined by accretion upon and radiation from the gainer. However, the radiation luminosity is much larger than the accretion luminosity. The luminosity of very large permanent disks exceeds the luminosity of the gainer at the end of RLOF. None of the initial systems has a period longer than 2.7 days and our disk systems all have observed periods between 4.5 and 13 days. This reflects the fact that the systems evolved in a conservative or weakly liberal way. Low values of the present mass ratio indicate that Algols that underwent RLOF during hydrogen core burning of the donor are at the end of their lives. Algols that underwent their first RLOF during hydrogen shell burning of the donor are in the middle of their RLOF lives despite the low value of their mass ratio.

Acknowledgments

We thank Monique and Jean-Louis Tassoul for their advice concerning tidal interaction and Norbert Langer giving us a copy of the Ph.D. thesis of Stephan Wellstein and for his suggestion to use a spectroscopic HRD in order to evaluate the mass dependent values of the overshooting parameter. We also thank Mercedes Richards, who is no longer with us, for applying Doppler tomography for the determination of disk characteristics.

References

All Tables

All Figures

Fig. 1 Evolution of SYNC = $\frac{{\mathit{\omega }}_{\mathrm{spin}}}{{\mathit{\omega }}_{\mathrm{orb}}}$, during the main sequence life of the 5 M⊙ primary in a 5 M⊙ + 2 M⊙ binary with a 10-day orbital period. From bottom to top: without tide, Darwin tide contributing almost nothing to the synchronization, weak Tassoul tide added, strong Tassoul tide added. The 5 M⊙ gainer needs the strong Tassoul tide to remain sufficiently synchronized.

In the text

	Fig. 2 Fractional radii (radius/semimajor axis) of six gainers and their accretion disks in the Lubow & Shu diagram. The observed radii of gainers are dots below ϖd. The calculated disk radii are squares above ϖd. The two observed disk radii are triangles. The fractional radii of permanent disks and the position of the gainer’s Roche radius are also indicated.
In the text

	Fig. 3 Evolution of the luminosity of the short-lived disk of β Lyr divided by the gainer’s luminosity as a function of relative disk time. From bottom to top: accretion luminosities (14), (18), radiative luminosities (24) and (20) added.
In the text

	Fig. 4 Evolution of the luminosity of the long-lived disk of TT Hya divided by the gainer’s luminosity as a function of relative disk time. From bottom to top: accretion luminosities (14), (18), radiative luminosities (24) and (20) added.
In the text

	Fig. 5 Evolution of the gainer in β Lyr through the HR-diagram. The short-lived disk is the thick part on the evolutionary path. The best fit model is the dot on the thick path whereas the observed position of the gainer is the dot to the lower left.
In the text

	Fig. 6 Evolution of the gainer of TT Hya. The dot on the thick path is the best model, whereas the observed position of the gainer is to the right.
In the text