Is the Blazhko Effect the Beating of a NearResonant DoubleMode Pulsation?
Abstract
In this paper it is shown that the Blazhko effect may result from a nearresonant type of multimode pulsation, where two (or sometimes more) periodic oscillations with slightly different frequencies gradually slip in phase, producing a beat frequency type of modulation. Typically one of these oscillations is strongly nonsinusoidal. Two oscillations are sufficient for the standard Blazhko effect; additional oscillations are needed to explain multifrequency modulation. Previous work on this hypothesis by Arthur N. Cox and others is extended in this paper by developing a simple (nonhydro) model that can accurately reproduce several important features found in Kepler data for RR Lyr, including the pulsation waveform, the upper and lower Blazhko envelope functions and the motion, disappearance and reappearance of the bump feature. The nonsinusoidal oscillation is probably generated by the fundamental mode and the other oscillations are probably generated by nonradial modes. This model provides an explanation for the strong asymmetry observed in the side peak spectra of most RR Lyrae stars. The motion and disappearance of the bump feature are shown to be an illusion, just an artifact of combining the oscillations. V445 Lyr is presented as an example with dual modulation. The mysterious doublemaxima waveform observed for this star is explained, providing additional strong evidence that this beatingmodes hypothesis is correct. Problems with other recent explanations of the Blazhko effect are discussed in some detail.
1 Introduction
Variable stars of the RR Lyrae type often exhibit a slow modulation of their pulsations known as the Blazhko effect. First reported in 1907 by Russian astronomer Sergey Nikolaevich Blazhko (Blazhko, 1907), the effect remains mysterious, with none of the numerous proposed explanations achieving widespread acceptance. See the recent reviews by Kovacs (2009) and Kolenberg (2012). Recent attempts to explain the Blazhko effect include:

“The shock model” in which a shock wave periodically generates turbulence in the atmosphere which in turn modulates the fundamental (Gillet, 2013).
All of these models have significant problems which are discussed in Section 4.
In this paper a simple explanation is put forth for the Blazhko effect which involves (at least) two periodic oscillations, PO1 and PO2, of slightly different frequency. Typically one of these oscillations is strongly nonsinusoidal. The slow phase slip between these oscillations results in a combined waveform that displays a beat frequency or apparent modulation at the Blazhko frequency. We refer to this as “the beatingmodes model”. It is the only model where there is no actual modulation of the amplitude of the primary mode, rather the apparent modulation is an illusion produced by the partial cancellation that results from adding to it the other mode (or modes) which have nearly the same frequency but differing phases. This hypothesis is supported by a simple (nonhydro) model of the dynamics that accurately reproduces Kepler data for RR Lyr, as shown in Figure 1. In addition, the beatingmodes model provides an explanation for the mysterious “doublemaxima” waveform of V445 Lyr (See Section 3.1), providing additional strong evidence in its favor.
Referring to Figure 1, note the many features of the data that are reproduced by the model including the shape of the pulsation waveform, the phase offset between the upper and lower modulation envelope functions, the size of the bump feature, the behavior of the bump: hesitating near the bottom of the pulsation waveform, moving upwards for a while, fading out and disappearing for a while and then reappearing near the bottom to complete the cycle. The disappearance of the bump lasts for about a third of the Blazhko cycle starting near the Blazhko minimum. The bump behavior of the model is a little surprising since it is entirely brought about by the addition of PO2, a pure sinusoid. PO1 by itself displays a bump at a fixed location about one quarter up from the bottom.
PO1 and PO2 may be generated by two unstable modes of the star. It is generally believed that the pulsation frequency for Blazhko stars corresponds to the fundamental mode, so this would be the likely choice for PO1, while PO2 may be generated by a nonradial mode of nearly the same frequency. An alternate possibility is that PO2 is a stable nonradial mode which becomes unstable by acquiring additional energy from the fundamental. (This is only possible if a mechanism exists that can accomplish this without phaselocking.) Since nonradial modes are involved, the usual 1D radial hydrocode software will not be able to test this model. MultiD codes are now being developed which may be able to accomplish this at some future date (see e.g. Mundprecht et al., 2013; Geroux & Deupree, 2014).
There have been several previous papers in support of the beatingmodes model, but the current work goes much further in terms of finding a simple model that actually provides a remarkably good fit to the observed Blazhko data, including the mysterious motion and disappearance of the socalled bump feature. Guggenberger et al. (2008), building on earlier work by Breger & Kolenberg (2006), proposed a test of the phasing behavior in an effort to determine if beating was responsible for the Blazhko effect. They applied this to the two wellobserved Blazhko RR Lyrae stars SS For and RR Lyr, with results consistent with the beating of two oscillations. Cox (1993), based on a combination of nonradial modal stability analysis combined with radial hydrocode analysis, proposed a nearresonant beating process between the fundamental and a nonradial mode, most likely the , mode. (Note that the stability analysis was for a static star and so might be significantly impacted by the presence of the high amplitude fundamental.) Sixteen years later, Cox (2009), working with improved code, finds that nonradial modes near the fundamental are always just slightly stable, and proposes that the , mode becomes unstable through nonlinear interaction with the fundamental. He still refers to this as a beating process not as a resonant mode interaction. Borkowski (1980) proposed a nearresonant doublemode pulsation involving the fundamental and either the second or third radial overtone, where the overtone frequency is higher than the second harmonic of the fundamental by the Blazhko frequency. Because the overtone is near the second harmonic instead of the fundamental itself, the Borkowski model is not validated by the results of this paper. Also since it only involves radial modes, one might have expected it to have been observed in 1D hydrocode models, but apparently it has not.
A variant of the beatingmodes model is found in the work of Nowakowski & Dziembowski (2001). Their approach involves a resonance between the fundamental mode and a pair of nonradial modes with opposite orders . Modes with nonzero have oscillations as a function of the azimuthal angle. These can be thought of as traveling wave modes that move around the star either with or against the rotational motion of the star depending on the sign of . As a result, these are split in frequency by the rotation. These modes can derive energy through nonlinear interaction with the fundamental if they satisfy the resonance condition that the sum of their frequencies must be exactly twice the fundamental frequency. As a result they generate a triplet spectrum with one peak on either side of the fundamental. According to Equation 77 of Nowakowski & Dziembowski (2001), the resulting Blazhko frequency is some fraction of the rotational frequency, most likely around 1/2. The pair are assumed initially stable and made unstable through the resonance. The resonant interaction between the modes is presumed stable and therefore the mode amplitudes and frequencies do not vary in time. Since all three peaks of the triplet correspond to actual modes, this is another type of beatingmodes model, where the combined waveform appears modulated even though no actual modulation taking place. The main problem with this model is that there seems to be no good explanation for the amplitude disparity that is often observed for the opposite side peaks. However it seems reasonable to consider, due to the wide variation of behavior observed among Blazhko stars, that this model could apply to some stars which have a relatively balanced triplet.
The spectra for Blazhko stars includes progressions of sidepeaks to the left and right of the fundamental and each of its harmonics, all evenly spaced apart by the Blazhko frequency. These tend to diminish in amplitude when moving away to the left or right and often only the first sidepeak is visible above the noise on each side. Most sidepeaks are nonlinear mixing product peaks and do not correspond to actual excited modes of the star. Usually these sidepeaks are attributed to modulation of the fundamental mode. However, if at least one of the sidepeaks immediately adjacent to the fundamental corresponds to an actual mode, this is sufficient through nonlinear processes (as discussed below) to generate the entire set of sidepeaks. The other frequencies are integer linear combinations of these two frequencies:
(1) 
where is a sidepeak frequency, and are integers (positive or negative), is the frequency of PO1 (the fundamental) and is the frequency of PO2. So, for example, the sidepeak that is opposite to (on the other side of the fundamental) has and . If one sidepeak is notably larger than the others then it is the likely choice for PO2. In the case of RR Lyr, there is a significantly larger peak on the higher frequency side. A spectral analysis (using the Period04 software package) of the Kepler data used in Figure 1(a) gives the amplitude (in units of ) of the upper peak as 0.0487 and of the lower peak as 0.0097, giving an amplitude ratio of about 5.02 and corresponding power ratio of 25.2. The central peak amplitude is about 0.1777.
It is well known that amplitude and frequency modulation can generate a sidepeak spectrum similar to the one observed (see e.g. Benko et al., 2011; Szeidl & Jurcsik, 2009; Szeidl et al., 2012), however a quasiperiodic dynamical system will also generate such a spectrum. A quasiperiodic system is one whose oscillations contain two base frequencies that are not rationally related. These would correspond to the frequencies of PO1 and PO2 and the apparent modulation frequency is the difference between them. The complete spectrum of such a system consists of the set of frequencies given by Equation 1. It can be obtained by a twodimensional Fourier series expansion as is shown in Appendix B. Nonlinear processes can change the amplitudes of the (sidepeak) frequencies including those which are initially zero. We take as a starting point a spectrum that includes PO1 (the fundamental mode frequency with its harmonics) and PO2 (a sine wave). There are two nonlinear mechanisms for generating the complete sidepeak array. One is the process that translates the pulsational motion of the star into a variation in the light output. There is no reason to assume this process is entirely linear. So, for example, if the light output depends slightly on the square of the sum of the radial velocities of the two oscillations, then there will be cross terms generated by the square which will introduce components at the triplet (or first order) sidepeaks for the fundamental and all of its harmonics. Higher order nonlinear terms will generate higher order sidepeaks. A second mechanism is generated by nonlinear coupling between the active modes. Although such coupling is (by assumption) insufficient to produce phaselocking between the modes, it can introduce sidepeaks that were previously missing, change the amplitudes of peaks that were already present and shift the base frequencies. In the first mechanism, the sidepeak frequencies are only present in the light curve, while in the second they are also physically present in the stellar vibration. Both mechanisms may be active at the same time.
2 Model
For purposes of matching the Kepler data, PO1 and PO2 are approximated as constant in amplitude and frequency throughout the Blazhko cycle, which is equivalent to saying that the nonlinear interactions between them are neglected as relatively small. PO1 is a highly nonsinusoidal waveform whose frequency is the pulsation frequency. Based on success at fitting the data with the model, it appears that this oscillation is made up of two subcomponents which we will call PO1a and PO1b both of which are sawtoothlike waveforms. PO1a is larger in amplitude, while PO1b is shifted in phase and generates the bump feature. The most likely possibility is that PO1 is generated by a single mode with complicated dynamics. It is generally thought that the bump is an echo, either a reflection from the stellar core or the atmosphere (Guggenberger & Kolenberg, 2006) and thus PO1b is the echo of PO1a. An alternate, but less likely, possibility is that PO1b is actually a separate mode, probably nonradial, that is phaselocked to the fundamental. This would make it similar to Cepheids which are thought to exhibit a bump due to a resonance with another mode (Simon, 1977; Buchler et al., 1990; Bono, Marconi & Stellingwerf, 2000). But regardless of the reason for the bump, the important point is that the model generates a reasonably good approximation to the observed waveform, enabling a test of the beatingmodes hypothesis. The model output for PO1 is shown in Figure 2. Details of how it is generated are given below.
PO2, on the other hand, is approximated as a pure sinusoid as this seems to do a good job of reproducing the observed Blazhko behavior. The frequency of PO2 is , where is the Blazhko frequency. In Section 5 it is discussed how the current work relates to the “hybrid mode model” developed previously by the author (Bryant, 2014).
For PO1a and PO1b we will use the same onezone model developed in Bryant (2014), although here we are only using the model to generate the sawtoothlike waveform needed for PO1. The model details are given in Appendix A. This onezone stellar model is similar to one developed by Baker (1966) and used extensively by Stellingwerf and others (Stellingwerf, 1972, 1986; Stellingwerf, et al., 1987; Munteanu et al., 2005) to study Cepheid and RR Lyrae variable stars. When excited to high amplitude, this model produces a sawtoothlike velocity function, a waveform that is commonly seen in the light curves of RR Lyrae and Cepheid variable stars. For simplicity, we will use identical model parameters to generate PO1a and PO1b, the only differences are that PO1b has its amplitude reduced by a scale factor and it is shifted in phase relative to PO1a. Both subcomponents are made equally nonsinusoidal, which seems to produce a very good fit to the observational data.
Let be the combined velocity of PO1 and PO2. For PO1 we use from Equation A4 in Appendix A and for PO2 we use a sine wave, with phase and scale factor chosen to fit the Kepler data:
(2) 
The light flux output , as plotted in Figure 1(b), is obtained from:
(3) 
The coefficients used were selected to obtain good visual agreement with the Kepler data in Figure 1(a)
We observe that for the Kepler data, the upper and lower Blazhko envelopes have differing amplitudes with the lower one being considerably smaller than the upper. The model on the other hand has upper and lower amplitudes of the envelope functions that are equal to each other and to the amplitude of PO2. We take a phenomenological approach to this problem and assume that the observed stellar flux is approximately equal to some nonlinear function of the combined radial velocity. In particular we add a correction term proportional to to the expression for the light flux. The result is shown in Figure 1(c) and obtained from
(4) 
The effect of the correction term is to stretch the upper part of the graph increasing the amplitude of the upper envelope and compress the lower part decreasing the amplitude of the lower envelope. Note that there may be other ways to achieve the same effect, e.g. it could conceivably result from a nonlinear interaction between the modes that generate PO1 and PO2.
3 Other Blazhko Cases
3.1 Multiple Frequency Blazhko
Recent work (Skarka, 2014; Benko et al., 2014) has shown that it is fairly common to find Blazhko cases that have two or more simultaneous modulation frequencies. This is accompanied with the appearance in the spectra of two or more sets of sidepeaks with different spacings. This could be the result of having an additional excited mode or modes near the fundamental that is not phaselocked to one of the peaks in the original sidepeak array. Nonlinear processes can then generate additional side peaks that are integer linear combinations of these with the other base frequencies.
A good example of a complex case is V445 Lyr (Guggenberger et al., 2012), which exhibits a dual frequency Blazhko modulation. The spectrum shows two strong peaks just above the fundamental and much weaker peaks below it. Thus it is likely that these strong peaks correspond to two active modes in addition to the fundamental. Each mode has an associated set of sidepeaks. Modeling in this case would require the addition of all three modal oscillations, which could be labeled PO1, PO2 and PO3. With this simple change, the model will generate the observed dual frequency Blazhko modulation. When the amplitude is large, the observed waveform has the typical sawtoothlike character (see e.g. Figure 1(b) in Guggenberger et al., 2012) although its rising section is not nearly as sharp as in the waveform for RR Lyr.
Occasionally the Blazhko minimum of V445 Lyr is very low, at which point it exhibits a “doublemaxima” waveform, see (see e.g. Figures 1(b) and 1(d) in Guggenberger et al., 2012). This can be easily explained as a cancellation of the base frequency, leaving the harmonic content behind. The first oscillation, PO1 is strongly nonsinusoidal and hence has a significant harmonic content. It can be thought of as the sum of two components: a sine wave at the base frequency and the harmonic content (2nd harmonic and above). The other two oscillations PO2 and PO3 have little harmonic content and so can be approximated as sine waves with frequencies that are very close to the base frequency of PO1. If the sum of the amplitudes of PO2 and PO3 is approximately equal to the amplitude of the base frequency component of PO1, then approximate cancellation is possible. It will occur when PO2 and PO3 are approximately in phase with each other and approximately 180 degrees out of phase with the base component of PO1. Since the relative phases of PO1, PO2 and PO3 are slowly changing, the near cancellation will only persist for a short time interval (a few pulsation cycles). The harmonic content of PO1 has nothing to cancel with and therefore remains, producing for that short interval a waveform that appears to be almost entirely harmonics. This is the observed doublemaxima waveform. The doubling of the maximum suggests a strong second harmonic content in PO1.
Explanations such as the ninth overtone model that involve an actual modulation of the fundamental would appear to be unable explain this effect. Normally when a nonlinear oscillation is at a low amplitude the nonlinear effects become dramatically less important. In the limit of small amplitude, the oscillation is expected to approach a perfect sine wave. But here the opposite is apparently happening—the amplitude is very small, but the harmonics, and particularly the second harmonic, are stronger than the base frequency.
Assuming for the sake of argument that the ninth overtone is a viable model (in spite of the serious problems outlined in Section 4.1), one might try to explain the strong harmonic content as resulting from interactions between the fundamental and the ninth overtone modes. But these interactions are nonlinear and in the small amplitude limit all such interactions vanish and the vibrational modes become completely independent. In this regard we note that the doublemaxima waveform appears only when the apparent amplitude of the fundamental mode is at it smallest, probably about 10 times smaller that at its maximum. We also note that the peak in the observed spectrum that (supposedly) corresponds to the ninth overtone mode is extremely small in amplitude, see in Table 2 in Guggenberger et al. (2012). Thus, within the ninth overtone model, it does not seem reasonable that sufficient nonlinear distortion can be present to account for the doublemaxima waveform.
3.2 Different Triplet Symmetry
A study of variable stars in the Large Magellanic Cloud (Alcock et al., 2003) found 74% of Blazhko stars had a larger sidepeak on the higher frequency side. These cases are potentially compatible with the RR Lyr results, but what about the other 26%? Since nothing presented here depends on the active sidepeak being on the high side, one possibility that the active mode is on the lower frequency side.
It is also possible that some stars may have additional sidepeaks that are active modes. If another mode is phaselocked to one of the sidepeaks (or to one of the harmonics of the fundamental) it may be difficult to determine with certainty that it is present. One indication of this could be that the sidepeak in question is notably higher in amplitude than expected from the pattern of its neighbors.
In cases where the side peaks are nearly symmetric in amplitude, the model of Nowakowski & Dziembowski (2001) could be at work. This is a beatingmodes model that is generated by a three mode resonance. See the discussion of this model in Section 1.
Cases of extremely low amplitude Blazhko modulation could potentially be generated by turbulence noise. This could be the case when there is a stable resonance between the fundamental and a nonradial mode at the same frequency. These modes would normally be locked at a particular relative phase. However, if the phase is offset from equilibrium it may oscillate back and forth as it decays back to the equilibrium point. Such oscillations are typically seen in amplitude equation simulations, see e.g. the introduction of Buchler & Goupil (1984). So there is phase modulation, but when adding the two modes together a beating process will occur generating the appearance of amplitude modulation. If driven constantly by noise, this might generate a symmetric triplet with fuzzy sidepeaks. This is of course rather speculative as it is not known if such phase locked states exist nor how much excitation could be achieved by the turbulence noise.
3.3 Long Period Blazhko
There are also cases with very long periods, up to 10 years (see e.g. Figure 5 in Soszynski et al., 2011), which seem to be at odds with the idea that nonlinear coupling should cause phaselocking for modes that are close together. One explanation is that the nonlinear coupling between the modes could be too weak for phase locking to occur. The coupling coefficient depends on the spacial and temporal characteristics of the modes and the order of the resonance involved and so it can vary substantially depending on the specific modes involved. There can also be an element of chance involved: If the coefficient is expected to lie in a range from negative to positive values then for certain stellar parameters it may “accidentally” have a value very near zero.
Another explanation is that these stars could have extremely low rotation rates. In this case, a dipolar nonradial mode could have its axis drift due to a combination of random forces from turbulence plus the effects of the slow rotation. If this mode was phaselocked to the fundamental, then this could be described as a noisy version of the Nowakowski & Dziembowski (2001) model (see discussion in Section 1) and the observed Blazhko frequency would be related to (but not equal to) the rotational frequency. As noted before, this model is also a beatingmodes model since the amplitudes of the modes are not actually being modulated.
This long period case is also problematic for Blazhko models of the “unstable resonance” type such as the ninth overtone model (Section 4.1). One might expect the system in that case simply go to stable phaselocking of the two modes rather than oscillate with an exceedingly long period—probably far longer than characteristic time constants of the system, such as the damping time constants for the modes involved. The time scale is also excessively long for the shock model (Section 4.2). It might be compatible with the Stothers model which seems to be more viable for very long period modulation (Section 4.3).
4 Problems With Other Blazhko Explanations
4.1 The Ninth Overtone Model
One explanation for the Blazhko effect (arguably the leading one) involves an unstable resonant interaction between the fundamental mode and the ninth radial overtone (Szabo et al., 2010; Buchler & Kollath, 2011; Kollath et al., 2011). This is a very high order resonance with the frequency ratio of the overtone to the fundamental being 9 to 2. As a result, correspondingly high powers appear in the amplitude equations (Equation 2 in Buchler & Kollath, 2011). This can cause extreme changes in coupling strength in response to moderate changes in amplitude.
Probably the most serious problem is that the spectral peak that they identify with the overtone is exceedingly small compared to the fundamental. Resonance based models like theirs depend on a nonlinear interaction between the modes to generate the Blazhko oscillation. This interaction must exhibit an instability that results in a substantial amount of energy being transferred back and forth at the Blazhko frequency between the fundamental and the other mode or modes. It is this energy shift that causes the observed amplitude changes in the fundamental during the Blazhko cycle. When the energy is shifted away from the fundamental it goes into the ninth overtone. One can see this quite clearly in the mode simulations of Buchler & Kollath (2011). In their Figure 1, note that the amplitude curve of mode B is roughly an inverted copy of mode A and is comparable in size. Also Buchler & Goupil (1984) state in their introduction “… one can have periodic or ‘chaotic’ energy transfer between the resonant modes.” But this appears to be inconsistent with the negligibly small size of the ninth overtone peak in the observational spectrum. To be specific, from Szabo et al. (2010) Table 1 and Figure 8, the amplitude of the peak at (the fundamental) is about 158 mmag and the peak at (the presumed ninth overtone) is about 1.65 mmag, about 100 times smaller. (The size discrepancy is even worse than it appears, since the kinetic energy in a mode is proportional to the square of the velocity, and the luminosity is known to vary roughly in proportion to velocity, making the energy ratio on the order of 10000.) When there is a large size discrepancy between the two modes there is a corresponding imbalance in the interaction between them. This is clear from the amplitude equations of Buchler & Kollath (2011) (their Equation 2) which shows that when is much smaller than , the effect of on the dynamics of is negligible, while the effect of on the dynamics of is quite strong. To make an analogy, an ocean wave can have a big effect on the motion of a surfer, while the surfer has a negligible effect on the wave. Thus the dynamics of the fundamental will proceed as if the overtone did not exist and there can be no Blazhko effect.
A second problem, is that the ninth overtone is being identified as a half integer harmonic of , but in addition to the expected peak at the spectrum also shows peaks at , , , , etc. Of these, peak at is the largest by far, and is close to 3 times the amplitude of the one at . This suggests that is the peak associated with an excited mode, while the others, including the peak at , may be nonlinear mixing product peaks. The fact that the peak is slightly elevated compared to the peak may simply reflect the fact that is the third harmonic of .
A third problem is that the alternating peak heights in the Kepler time series data appear and disappear in an intermittent, seemingly random fashion, rather than varying in synchronization with the Blazhko cycle. If the half integer harmonics are directly related to the Blazhko effect then they must play an active role in the Blazhko dynamics. They cannot be exhibiting their own separate dynamics and still be held responsible for the Blazhko effect. This is in contrast to the Smolec et al. (2011) test of the Stothers model, where they state “We note that in our models period doubling is strictly repetitive and always appears at the same phase of the Blazhko cycle.”
A fourth problem is that since the ninth overtone model is entirely radial, the Blazhko effect should be reproducible with the available 1D hydrocodes. Apparently this has not occurred.
4.2 The Shock Model
The shock model of Gillet (2013) presents just the framework of a possible theory with all of the needed details either missing or fatally flawed. The basic idea is that shock waves build up in the atmosphere while the fundamental is growing in amplitude, then these shocks generate turbulence affecting the excitation mechanism for the fundamental, then the fundamental decays in amplitude and finally the atmosphere becomes stable again and the process repeats. Gillet claims this process depends critically on the presence of the “transient” first overtone, but this concept is not clearly defined. Does transient mean the mode is only present for part of the Blazhko cycle? If so what portion of the cycle and why does Gillet not look at the observational data to try and prove that the overtone has a transient presence? Gillet admits that the observed amplitude of the first overtone is about 100 times smaller than the fundamental. From this one might estimate that the energy could be on the order of 10000 times smaller since energy typically is proportional to the square of the velocity while the light signal tends to be roughly proportional to the velocity. So how can something so weak and insignificant be expected to play a major role in the stellar dynamics? Gillet also admits that “The physical origin of shock s3’ has not been determined yet” (i.e. it is complete speculation that is has anything to do with the first overtone). Also there is the fact that the first overtone has a different frequency than the fundamental so if it is supposed to be generating a shock wave at a certain phase of the fundamental pulsation, there will be a problem since the overtone will have a different relative phase each successive pulsation. The results also depend on the ability of the shock wave to suddenly cause the complete “desynchronization” of the photospheric layer; is there evidence that a shock wave can in fact do this and do it suddenly rather than gradually? One should bear in mind that these shock waves are not caused by an outside source, but simply appear as part of the dynamics of the pulsational motion. So it seems surprising that such a shock wave would suddenly disrupt the motion that is generating it. It would seem more likely that if turbulence was emerging this would limit the strength of the shock wave and an equilibrium would be reached. Aside from the turbulence associated with the desynchronization this would seem to be an entirely radial model and so there is also the question as to why the proposed mechanism is not seen in 1D hydrocode simulations. Fossati et al. (2014) discuss this model—they find that the changes in the level of turbulence required by the model appears higher than what should be expected from their measurements of variations in the microturbulent velocity for RR Lyr.
4.3 The Stothers Model
The Stothers model has received a fair amount of attention recently. It is discussed at length in Kovacs (2009) and the papers by Smolec et al. (2011) and Molnar et al. (2012) are entirely devoted to it. It depends on several assumptions that are very difficult to test. In particular it is assumed that there is a magnetic field whose amplitude varies at the Blazhko frequency and that the changes in magnetic field in turn affect the level of turbulence sufficiently to modulate the amplitude of the fundamental mode. In an effort to test the theory, Smolec et al. (2011) use a hydro model into which the assumed periodic changes in the turbulence are introduced in an ad hoc manner. They are able to recreate some of the desired Blazhko features including asymmetric sidepeaks and the disappearance of the bump during the ascending Blazhko phase, but they think that the amount of turbulence changes required for this may unrealistic. The simulation also fails to reproduce the upward motion of the bump as seen in Figure 1 of the current paper. The behavior of the bump in the Stothers case is of course a function of the level of turbulence. Also there is no significant amplitude differential between the upper and lower Blazhko envelope functions. Molnar et al. (2012) using a different hydro model find the Stothers model could be workable but only for Blazhko periods longer than 100 days due to the slow response of the system. Fossati et al. (2014) discuss this model—they find that the changes in the level of turbulence required by the model appears higher than what should be expected from their measurements of variations in the microturbulent velocity for RR Lyr.
5 Relation to Earlier Work by the Author
In a previous paper (Bryant, 2014) it was shown that an accurate fit to second quarter Kepler data for RR Lyr could be achieved by a hybrid mode consisting of two component modes of the same frequency. The first of these components is a highly nonsinusoidal sawtoothlike waveform, while the second is approximated as a sinusoid. Although locked in frequency, the two components are allowed to change independently as a function of the Blazhko phase. The first component varies little in amplitude or frequency throughout the Blazhko cycle, while the second changes quite strongly both in amplitude and in phase relative to the first component. By optimizing the fit to the Kepler data, the amplitude and phase behavior was obtained and displayed in Figure 4 of Bryant (2014).
The current work is related to that earlier work in the following way: It is found that the phase and amplitude dynamics of the second component of the earlier description can be approximately generated by adding together two constant amplitude sine waves, one with frequency and a second of lower amplitude with frequency . This was noted previously by Guggenberger et al. (2008) who found this type of phase and amplitude behavior in several RR Lyrae stars and attributed it to this type of beating process. The first component corresponds to PO1a in the current model. The second sine wave corresponds to PO2. The first sine wave corresponds to PO1b, although in the current model this is made nonsinusoidal. Making this nonsinusoidal helps to correct one of the problems with that earlier model, in that the bump feature generated by that model was not as pronounced as in the Kepler data.
6 Conclusion
It has been shown that a simple beatingmodes model of the Blazhko effect can provide a remarkably good fit to Kepler data for RR Lyr. The model consists of two periodic oscillations of slightly different frequency that are simply added together. This is in contrast to the much more complex picture of resonant mode models that depend on nonlinear coupling between modes leading to an oscillatory exchange of energy between them. The current model reproduces for the first time the apparent motion, disappearance and reappearance of the bump feature, the phase difference between the upper and lower modulation envelope functions and their relationship to the motion of the bump. In addition, it has been shown that the amplitude disparity between the upper and lower envelope functions can be reproduced if the light curve is assumed to be a nonlinear function of the radial velocity. Previous work by Cox (1993, 2009) suggests that the modes involved could be the fundamental and the , nonradial mode. The model is easily extended to cases exhibiting multiple modulation frequencies by simply associating additional large nearresonant peaks in the spectrum with active modes of the star, i.e. one additional mode is required for each additional modulation frequency. Application of this extended model was discussed for the case of V445 Lyr. The mysterious doublemaxima waveform observed for that case is explained providing very strong evidence that the hypothesis of this paper is correct.
In addition, problems with other recent models of the Blazhko effect have been discussed in some detail, including the ninth overtone model, the shock model and the Stothers model.
Appendix A OneZone Model and Po1 Details
The equations of motion are derived from a simplified analysis of pressure and gravity acting on the portion of a star above a certain radius . We may take to be the effective bottom of the stellar envelope if we are considering the fundamental mode (as we are in this paper) or the outermost radial node for other modes such as the first radial overtone. For large amplitude oscillations, the motion is ballistic or gravity dominant for most of the cycle, with pressure becoming dominant only near the minimum position of the cycle. This leads to a “bouncing ball” type of oscillation whose derivative is the sawtoothlike waveform mentioned previously. Regardless of the accuracy of this simplified picture of the dynamics, the ability to reproduce the sawtoothlike waveform shows that the onezone model seems to capture an essential aspect of the true dynamics.
The model takes the form of an ordinary differential equation for a harmonic oscillator, but with the usual linear forcing term replaced with a nonlinear function that, in a very simplified way, captures the essence of the gravity and pressure forces, i.e. and where and are position and velocity. This model is conservative, excitation and damping of the oscillations are omitted, but could be added to the model if desired. Also currently assumed to be small and ignored are nonlinear interactions between modes that could lead to energy transfer or phase shift. In Bryant (2014) the following expression is derived for the forcing function:
(A1) 
where:
(A2) 
(A3) 
, , is a characteristic radius for the stellar material that lies above and is the polytropic index. The problem has been rescaled so that . To keep the values near unity, time is in units of 0.1 days. The values used in the model are rad/0.1 day, , , and . The Blazhko period is 72 pulsation periods. Initial condition for PO1a: , . Initial condition for PO1b is chosen so that it produces the identical oscillation except that it is leading in phase by 2.1058232 radians or 0.19 days: , . Note that is the frequency for small amplitude oscillations; the actual oscillation frequency drops as the amplitude increases and is equal to for the specified initial conditions. The velocity for PO1, , is the combined velocity of PO1a and PO1b, with scale factors chosen to fit the Kepler data:
(A4) 
This function is plotted in Figure 2.
Appendix B Frequency Spectra for Quasiperiodicity
When a nonlinear dynamical system is quasiperiodic, it will have two base frequencies which are incommensurate, i.e. they do not have a rational ratio. In the phase space the dynamics will lie on a twotorus or doughnut surface. The dynamics will wind its way around this surface, taking an infinite amount of time to completely cover it. Motion on a 2torus can be described by two anglelike variables, and , corresponding to the two ways in which one can go around the torus. We can define the angles so that they satisfy and . The state of the system can be expressed as a function of these two angles:
(B1) 
The state is periodic in both angles (since a shift in either angle corresponds to one pass around the torus). Thus it can be expanded as a two dimensional Fourier series:
(B2) 
where
(B3) 
In terms of time this becomes:
(B4) 
where the frequencies in the spectrum are
(B5) 
The amplitudes can also be obtained as an integral over time:
(B6) 
The amplitude of these components must fall off as and become large.
Cases of three frequency quasiperiodicity (or higher) require a simple extension of the above analysis to a higher dimensional Fourier series.
References
 Alcock et al. (2003) Alcock, C., Alves, D. R., Becker A., et al. 2003, ApJ., 598, 597
 Baker (1966) Baker, N. 1966, in Stellar Evolution: Proc. of International Conf., November 1315, 1963, ed. R. F. Stein & A. G. W. Cameron eds., (New York: Plenum), 333
 Benko et al. (2011) Benko, J. M., Szabo, R. & Paparo, M. 2011, MNRAS 417, 974
 Benko et al. (2014) Benko, J. M., Plachy, E., Szabo, R., et al. 2014, ApJ Supplement Series, 213, 31
 Blazhko (1907) Blazhko, S. N. 1907, Astron. Nachr., 175, 325
 Bono, Marconi & Stellingwerf (2000) Bono, G., Marconi, M. & Stellingwerf, R. F. 2008, A&A 2000, 360, 245
 Borkowski (1980) Borkowski, K. J. 1980, Space Science Reviews 27, 511
 Breger & Kolenberg (2006) Breger, M. & Kolenberg, K. 2006, A&A 460, 167
 Bryant (2014) Bryant, P. H. 2014, ApJ 783, L15
 Buchler & Goupil (1984) Buchler, J. R. & Goupil, K. 1984, ApJ 279,394
 Buchler & Kollath (2011) Buchler, J. R. & Kollath, Z. 2011, ApJ 731, 24
 Buchler et al. (1990) Buchler, J. R., Moskalik, P. & Kovacs, G. 1990, ApJ, 1990, 351, 617
 Cox (1993) Cox, A. N. 1993, in Proc. IAU Coll. 139, New Perspectives on Stellar Pulsation and Pulsating Variable Stars, Eds. J. Nemec , J.M. Matthews, Cambridge Univ. Press, p. 409
 Cox (2009) Cox A. N. 2009, in AIP Conf. Ser. Vol. 1170, Stellar Pulsation: Challenges for Theory and Observation. Am. Inst. Phys., New York, p. 276
 Fossati et al. (2014) Fossati E., Kolenberg, K., Shulyak, D. V., et al. 2014, MNRAS Accepted, arXiv:1410.4380 [astroph.SR]
 Geroux & Deupree (2014) Geroux, C. M. & Deupree, R. G. 2014, ApJ, 783, 107
 Gillet (2013) Gillet, D. 2012, 2013, A&A, 554, 46
 Guggenberger & Kolenberg (2006) Guggenberger, E. & Kolenberg, K. 2006, Comm. in Asteroseismology., 148, 21
 Guggenberger et al. (2008) Guggenberger, E., Kolenberg, K. & Breger M. 2008, Journal of Physics: Conference Series, 118, 012056
 Guggenberger et al. (2012) Guggenberger, E., Kolenberg, K., Nemec, J. M., et al. 2012, MNRAS 424, 649
 Kolenberg (2012) Kolenberg, K. 2012, JAAVSO, 40, 481
 Kolenberg et al. (2010) Kolenberg, K., Szabo, R., Kurtz, D. W., et al. 2010, Astrophys. J. Lett., 713, L198
 Kollath et al. (2011) Kollath, Z., Molnar, L. & Szabo R. 2011, MNRAS 414, 1111
 Kovacs (2009) Kovacs, G. 2009, AIP Conference Proceedings 1170, 261
 Molnar et al. (2012) Molnar, L., Kollath, Z. & Szabo R. 2012a, MNRAS, 424, 31
 Mundprecht et al. (2013) Mundprecht, E., Muthsam, H. J. & Kupka, F. 2013, MNRAS, 435, 3191
 Munteanu et al. (2005) Munteanu, A., Bono, G., Jose, J., et al. 2005, ApJ, 627, 454
 Nowakowski & Dziembowski (2001) Nowakowski, R. M. & Dziembowski, W. A. 2001, Acta Astronomica, 51,5
 Simon (1977) Simon, N. 1977, ApJ, 1977, 217, 160
 Skarka (2014) Skarka, M. 2014, A&A 562, A90
 Smolec et al. (2011) Smolec, R., Moskalik, P., Kolenberg, K., et al. 2011, MNRAS, 414, 2950
 Soszynski et al. (2011) Soszynski, I., Dziembowski, W. A., Udalski, A., et al. 2011, Acta Astron., 61, 1
 Stellingwerf (1972) Stellingwerf, R. F. 1972, A&A, 21, 91
 Stellingwerf (1986) Stellingwerf, R. F. 1986, ApJ 303, 119
 Stellingwerf, et al. (1987) Stellingwerf, R. F., Gautschi, A., & Dickens, R. J. 1987, ApJ, 313, L75
 Stothers (2006) Stothers, R. B. 2006, MNRAS 365, 178.
 Stothers (2010) Stothers, R. B. 2010, PASP 122, 536.
 Szabo et al. (2010) Szabo R., Kollath, Z., Molnar, L., et al. 2010, MNRAS, 409, 1244
 Szeidl & Jurcsik (2009) Szeidl, B. & Jurcsik, J. 2009, Comm. in Asteroseismology, 160, 17
 Szeidl et al. (2012) Szeidl, B., Jurcsik, J., Sodor, A., Hajdu, G. & Smitola, P. 2012, MNRAS 424, 3094