On the Modeling of Algol-Type Binaries

: In earlier papers, we presented a binary evolutionary code for the purpose of reproducing the orbital parameters, masses, radii, and location in the Hertzsprung Russell diagram (abbreviated as HRD) of well-observed Algol systems. In subsequent versions, the effects of mass and angular momentum losses and tidal coupling were included in order to produce the observed distributions of orbital periods and mass ratios of Algol-type binaries. The mass loss includes stellar wind and possible liberal evolution, when the gainer star is not capable to absorb all of the matter during mass transfer from the donor star. We added magnetic braking to our code to better reproduce the observed equatorial velocities. Large equatorial velocities of mass-gaining stars are now lowered by tidal interaction and magnetic braking. Tides are mainly at work at short orbital periods, leaving magnetic braking alone at work during longer orbital periods. The observed values of the equatorial velocities of mass gainers in Algol-type binaries are mostly well reproduced by our code. According to our models, Algols have short periods with a strong magnetic ﬁeld. rotational velocities (close to the critical one) and the observed values. In this paper, we present a brief overview of the improved modelling work that has been done over more than a decade, addressing these issues with some updated results.


Introduction
The modelling of the evolution of close binaries started with the first papers of a series from Paczynski [1][2][3], soon followed by papers from colleagues in Germany (Kippenhahn & Weigert [4]) and from the USSR (Tutukov & Yungelson [5,6]). Paper [5] explored the loss of matter and angular momentum during the period of mass transfer (further called Roche Lobe Overflow, which is abbreviated as RLOF). Further developments can be found in the review that was published by Eggleton et al. [7]. Nelson & Eggleton [8] used a fast code producing large amounts of models to allow statistical studies. However, the code was less suitable to produce models matching individual systems due to the use of averaging equations. Song & Maeder et al. [9] and Song & Meynet et al. [10] provided fundamental refinements of the theory of evolution of interacting binaries.
Van Rensbergen & De Greve et al. [11] started a series of papers for explaining the observed distributions of mass ratios and orbital periods of Algol-type systems. Gradually, improvements were introduced in the code to obtain better matches between observations of individual systems and models. In the papers [12,13], we investigated mass loss from the system (the liberal case) during eras of rapid mass transfer, while using a hot spot mechanism. In [14], we presented models for systems with an accretion disk. Equatorial velocities of gainers were modelled in [15] using magnetic braking to keep the equatorial velocity far below critical during the evolution of the binary. These successive studies still revealed some remaining discrepancies between theoretical and observed mass ratio and orbital period distributions and raised new questions on the difference between the theoretical rotational velocities (close to the critical one) and the observed values. In this paper, we present a brief overview of the improved modelling work that has been done over more than a decade, addressing these issues with some updated results.

Paths through the Hertzsprung Russell Diagram (HRD)
Every observed binary has a large number of possible progenitors with different values for their masses and different initial periods. Our evolution code recognizes a progenitor as plausible when the calculated positions of donor and gainer in the HRD are found to be close to the observations. The HRD is an effective temperature-luminosity diagram. Observers refer to it as colour -magnitude diagram (abbreviated as CMD). As an example, Figure 1 shows the evolutionary path of β Per: prototype of the Algol-type binaries. The present characteristics of Algol are from the catalogue of Budding et al. [16]. All characteristics have been accurately retrieved by our code, except the T eff -values. However, Zavala et al. [17] revised these T eff -values. Their value of T eff for the gainer is closely approximated by our calculations. The progenitor is a (3.41 M + 1.1 M ) binary with an initial period of 1.4625 days. This binary lives two eras of RLOF. Nowadays, the second era of RLOF (RLOF B during H-shell burning of the donor) is at work. Figure 1 shows that the observed location of the gainer is produced by the calculations, whereas the calculated effective temperature of the donor remains larger than the observed one, as mentioned in the catalogue of Budding et al. [16]. Figure 1 also shows the locations on the Zero Age Main Sequence, where both stars start their life. During evolution, the originally most massive donor becomes the less massive member of the binary. Presently, the most massive gainer is still on the Main Sequence (abbreviated as MS) of the HRD, whereas the less massive donor has evolved far away from the MS. The Algol paradox states that the most massive member of an Algol-binary is less evolved than his less massive companion. Figure 1 shows the well-known solution to this paradox.

Tidal Interaction
Neglecting tides and magnetic braking, Packet [18] showed that mass gainers continue to rotate with the critical velocity after the accretion of only 5% to 10% of their own mass. This critical velocity is, in most cases, much faster than that observed. Van Rensbergen et al. [11] calculated a grid of evolutionary tracks with a simple approach for tidal interaction and the possibility for binaries to live through liberal eras during their evolution.
Darwin [19] developed tidal interaction as first. Tides act on both stars of the binary, depending on their radius "R" and the semi major axis "a" of the binary orbit. A typical synchronization time is given by: The mass of the star undergoing the tide is in the denominator of the mass ratio q in Equation (1). Tides tend to synchronize the angular velocities of the stellar rotation (Ω) with the angular velocity of the orbit (Ω orb ), as shown in Equation (4).
In order to distinguish strong tides from weak tides, Wellstein [20] proposed calculating the synchronization time with: t sync = τ sync f sync (2) For strong tides, Wellstein [20] uses f sync = 0.1, for weak tides f sync = 1 is taken in Equation (2). Both of the values of f sync were used by Van Rensbergen & De Greve et al. [11][12][13]. Strong tides should be preferred over weak tides in order to obtain a better agreement between observed and calculated distributions of mass ratios an orbital periods. Nevertheless, strong tides are also not able to predict the population of Algol-binaries with the smallest observed orbital periods.
In order to avoid the random number f sync , Van Rensbergen & De Greve [14] calculated t sync with a physical model, distinguishing convective from radiative envelopes. Hurley et al. [21] use an expression of Hut [22] to calculate t sync for a star with a convective envelope. Hilditch [23] gives an expression of Zahn [24] to calculate the tidal action on a star with radiative envelope. In the following, we will designate these synchronization times with t sync,Darwin . Assuming that mass accretion causes turbulence in the atmosphere, the convective mode was always applied during RLOF. Meridional circulation, as proposed by Tassoul [25], was added as a significant contributor to tidal interaction. Figure 2 shows the effect of Tassoul's mechanism, presenting the degree of synchronism of the 8 M star in a (8 M + 3 M )-system in the pre-RLOF phase. During RLOF, the separation of different tidal effects is blurred by the mass transfer. Thus, the tidal synchronization time was calculated with: Tidal synchronization tries to equalize the angular velocities of the stars (Ω) in a binary with angular velocity (Ω orb ) of orbital rotation. Because of tides, a rotating star with moment of inertia I and angular velocity Ω will change its spin angular momentum by an amount This equation shows that tides cease their action on the rotation of a star when synchronism is achieved with Ω = Ω orb .

Liberal Evolution
Peters & Polidan [26] introduced the concept of hot spots (also called HTAR: High Temperature Accretion Regions) that are created at the trailing side of the gainer by the mass-transfer from the donor hitting the gainer star. Van Rensbergen & De Greve [12,13] used these hot spots to act together with rapid rotation of the gainer, attempting to surpass the binding energy at the HTAR-location. More massive binaries reach that condition, as they yield both characteristics, fast rotation and a very hot spot, during periods of rapid mass transfer. The evolution of the binary will be liberal during such an era. The influence of liberal evolution on the distribution of orbital periods and mass ratios of the evolutionary models was compared to the observed distributions for 351 Algols, which were taken from the catalogue of Budding et al. [16]. This catalogue uses different values of the mass ratio: q LC is the mass ratio that is derived with the light curve solution and q MS is the mass ratio derived to make the parameters of the gainer stars to fit Main Sequence characteristics. Figure 3 compares observations and our models for the mass ratio distribution q = M d M g , where the subscript "d" is used for the donor star and "g" for the gainer. We updated the comparison by taking the quantity q LC for the observed distribution instead of q MS , as the latter value results in larger uncertainties than the former. The agreement between observations and model is better than satisfactory, as compared to the result that was obtained by Van Rensbergen et al. [13], with the injudicious use of q MS , which produces too many Algols with large mass ratios, a result that is contradicted by observations. The majority of Algols, having lower mass progenitors, evolve conservatively. The less numerous more massive Algols evolve very differently as a consequence of their liberal eras. However, since most Algols evolve conservatively, conservative and liberal simulations yield similar distributions of the mass ratio.
The initial conditions for the simulations are from Salpeter [27] for the Initial Mass Function (IMF) of the future donor, from Van Rensbergen et al. [28] for the initial mass ratio distribution, and from Popova et al. [29] for the distribution of the initial orbital periods. Figure 4 shows that the observed distribution of orbital periods of 376 Algols from the catalogue of Budding et al. [16] is well reproduced by our simulation. The difference of 376 with the earlier mentioned 351 comes from the fact that the period is known for 25 systems, for which the mass ratio still remains to be determined. In this case there is no difference between the results obtained from the conservative and the liberal simulation. However, our simulations yield almost no Algols with orbital periods shorter than one day, in contrast with the ≈10 % found in the catalogue.

Magnetic Braking
Tides and liberal evolution cannot explain the observed rotational velocities of the mass gainers published by Van Hamme and Wilson [30], Miller et al. [31], Glazunova et al. [32], and Dervisoglu et al. [33]. Therefore, we investigated another possible process. In Section 5.1 we explain how to calculate the magnetic field strength, while using a simple model of differential rotation. The resulting magnetic braking rises with increasing values of the magnetic field and the amount of mass loss through stellar wind. The mass losses by stellar wind are from Vink et al. [34] for stars hotter than 12,500 K and De Jager et al. [35] for cooler stars. Details of the calculation method can be found in Van Rensbergen & De Greve [15]. We further report the results, including magnetic braking, in order to better understand the observed rotational velocities of Algol-type systems.

Generation of the Magnetic Field of the Gainer
A solid body rotator with constant angular velocity Ω from the centre to the edge of the star cannot develop a magnetic field. We use a model, wherein only the outer shell is spun up by the mass-transfer from the donor hitting the gainer. During rapid RLOF the outer shell rotates faster than the core of the gainer whose rotation is not accelerated. The magnetic fields that are produced by the dynamo model of Spruit [36] are at work when the angular velocity increases with an amount ∆Ω over a distance ∆r. Packet [18] only gives the amount of angular momentum that is added to the shell, being reduced by the quotient of the impact-parameter d over the radius of the gainer R g (rendering the spin-up inefficient at close approach): ∆J + spin,shell = 6.04534 10 54 R g M g + ∆J + spin,shell is expressed in cgs units, whereas masses and radii are, respectively, in M and R .
With every value of J spin,shell corresponds a value of Ω shell = J spin,shell I shell , which is a characteristic value of the angular velocity in the shell. A magnetic field is created when Ω shell > Ω core . However, the value of Ω cannot change discontinuously from Ω core to Ω shell at the interface between core and shell. In order to accomplish the continuity of Ω(r) and to avoid friction between the outer edge of the core and the inner surface of the shell, we chose a radius-dependent angular velocity Ω(r) continuously increasing from Ω core at the interface between core and shell to Ω edge at the edge of the gainer.

Extent of the Shell
The gainer is now composed of an inner core in solid rotation that was surrounded by a shell that rotates differentially. The core is not spun up or slowed down magnetically. Its rotation is only modulated by tides. When the shell is spun up following relation (5), the spin-down follows relation (4).
A core with a core fraction CF of the mass of the gainer leaves (1-CF) for the mass for the shell. Figure 5 shows the evolution with time of the equatorial velocity of the gainer starting from a (6.36 + 2.7) M binary with an initial orbital period of 1.89745 d. This is a plausible progenitor for λ Tau, although any other binary in our collection of calculated cases could have been chosen to illustrate our updated binary evolutionary scenario. The calculation was performed for the system undergoing tides and magnetic braking with an initial value of CF = 0.95. The evolution with time is compared with the evolution of the gainer in the same system with a rigidly rotating and, hence, never magnetic gainer, for which the rotation is only affected by tides. Also included in the figure is the evolution when neither tides nor magnetic braking are active, showing that, in the latter case, the gainer remains rotating with critical velocity once this velocity is achieved. The plausible progenitor for λ Tau starts RLOF during core hydrogen burning of the donor (RLOF A). A critical rotation of the gainer is achieved during this phase of rapid RLOF at ≈48 million years after ZAMS. After that, RLOF A occurs at a lower speed. Tidal interaction and magnetic braking are then strong enough to synchronize the rotation of the gainer. In Figure 5, one clearly sees that synchronization settles more rapidly when magnetic braking helps tidal interaction. After ≈65 million years, during hydrogen shell burning of the donor, RLOF restarts. Critical rotation is again attained. However, the rush to critical rotation is slowed down by the combined action of tides and magnetic braking. The present observed state of λ Tau occurs during this stage of rising rotational velocity. The gainer of λ Tau now has an equatorial velocity that is far below the critical value.
When critical rotation is achieved during RLOF B, the orbital period is ≈30 days. Figure 5 shows that tides are then too weak to prevent critical rotation. However, magnetic braking will decrease the velocity below the critical value. This happens at the very end of, and after, RLOF B. Figure 6 shows the evolution of the magnetic field of the gainer in the binary that is shown in Figure 5. Before RLOF A, the shell and core of the gainer have the same synchronous rotational velocity. The core always continues to rotate almost synchronously. When RLOF A starts, only the shell is spun up. When the shell rotates at critical velocity, the magnetic field is at a maximum of ≈3000 Gauss for a short time. The magnetic field disappears when tides and magnetic braking have synchronized the rotation of the shell. From the beginning of RLOF B, the shell is again spun up, and the magnetic field starts to be built up again. The value of the magnetic field is lower than during rapid RLOF A. A value of ≈750 Gauss is reached. The magnetic braking is then still sufficiently active to synchronize the shell. After that, the magnetic field disappears again.

Calculated Cases
Van Hamme & Wilson [30] defined the quantity F= v eq v sync : the ratio between the equatorial velocity v eq and the synchronous velocity v sync . Subsequently, they determined R= F gainer −1 F crit,gainer −1 . This quantity R ∈ [0-1] is a measure of rotation. Specifically, R is 0 at synchronous rotation and 1 at critical rotation. The use of the quantity R as a characteristic measure of rotation has the advantage that minor differences in the evaluation of the synchronous and critical velocities by different authors are ironed out. A determination of v eq is considered to be good when ∆R = |R model − R obs | is small (near to 0) and it is unacceptable when the same quantity is large (close to 1). We found 41 Algol-type binaries for which the equatorial velocity of the gainer was published. These systems are shown in Tables 1 and 2, running from ∆R = 0 at the top of Table 1 to the bottom of Table 2 where ∆R ≈ 1. The Tables show two lines for each system (first line: observed data; second line: data from the model), showing that: Only one Algol (U CrB) is badly reproduced by the theory (∆R > 0.90). Table 1. Calculated equatorial velocities of gainers that best fit observations. Masses are in M , orbital periods in days and v eq in km sec . For each system, two lines are given, the first with data from observations and the second with the results from the model. The masses mentioned by the observations of Glazunova et al. [32]

Conclusions
For a large part of Algol-type binaries with well observed characteristics, fairly accurate models can be calculated matching the present position in the HRD. For some systems, a better determination of fundamental parameters, such as T eff for the cooler donor star is needed. The observed distributions of mass ratios and orbital periods of Algol-type binaries are reasonably well reproduced by our models, although additional work is needed to represent the very short period Algols, as can be seen in Figures 3 and 4. We suggest shedding new light on the contact phase evolution, where very short orbital periods are the rule. The present calculations show that magnetic braking added to tidal interaction clearly helps to reduce the rotational equatorial velocity of the gainer star, away from the critical one. Reasonable to good agreement between the observed and calculated rotational velocities was obtained for 33 out of 41 gainers in Algol-type binaries. The unsatisfying result for eight systems requires further investigation. It may be due to different causes. The current period is well known, but that is not the case for the masses. If future analysis would deliver other values for the observed masses, then this would result in different progenitors, subsequently with a different v eq -evolution. Additionally, the eight systems all occur long after the minimum period (late case A or case B). A new approach towards their previous contact phase may also produce different results. Our models predict that Algol-type binaries may show a strong magnetic field for a short time during phases of rapid mass transfer.
Author Contributions: The contributions of two authors are equal. Both authors have read and agreed to the published version of the manuscript.