EFFECTS OF NON-CIRCULAR MOTIONS ON AZIMUTHAL COLOR GRADIENTS

Until the detection of an azimuthal color gradient across one of the arms of the SAc galaxy M 99 (González & Graham 1996, GG96 hereafter), only sparse evidence of star formation triggered by spiral density waves had been found in stellar counts in the Milky Way (Sitnik 1989, 1991; Avedisova 1989) and M 31 (Efremov 1980a, 1980b, 1985). More recently, applying the same method defined by GG96 and described below, Martínez-García et al. (2009, MG09 hereafter) examined a sample of 13 spiral galaxies of types A and AB, and found color gradients consistent with theoretical expectations in 10 of their objects.³ Although they did not compute stellar orbits, both GG96 and MG09 implicitly assumed that all stars, including those recently born near spiral arms, move in circular orbits. MG09 did investigate the effects of variable circular speeds and variable densities on observed color gradients and found them to be negligible. However, complicated non-circular motions have been reported in studies about the migration of young stars following star formation triggered by spiral shocks (Yuan 1969; Wielen 1979; Fernández et al. 2008). Bash (1981) had also noticed that the ballistic orbits calculated from galactic H ii region complexes follow non-circular trajectories that initially move along (and not across) the spiral arms.

In the present work, we examine how the presence of non-circular motions would modify the pattern speeds derived from color gradients under the assumption of circular orbits.

We will explore two different approaches. The first one (see Section 3) involves the semianalytical solutions obtained for spiral shocks (Roberts 1969; Shu et al. 1973; Gittins & Clarke 2004); we assume that newly born stars preserve the orbital motion of the shocked gas where they form. The second approach (see Section 4) is based on data from magnetohydrodynamic (MHD) simulations; we follow the gas flow vectors near spiral arms. We always assume that stars are triggered all along the studied regions of spiral arms in the respective model. Other methods, not discussed here, may involve orbit calculations for young stellar groups.

The photometric technique employed by GG96 and MG09 uses three optical bands, g, r, and i, plus the near-infrared J band (see Table 1). With these filters, the reddening insensitive and star formation sensitive photometric index Q is obtained:

$$\begin{equation} Q(rJgi) = (r-J) - 0.82 (g-i), \end{equation}\vspace*{-11pt} \tag{ 1 }$$

$$\begin{equation} Q(rJgi) = \log _{10} \frac{I^{2.05}_{g} I^{2.50}_{J}}{I^{2.50}_{r} I^{2.05}_{i}}. \end{equation} \tag{ 2 }$$

In star-forming regions, the Q(rJgi) index has higher values, because the g and J bands in the numerator of Equation (2) trace the light from red and blue supergiants, respectively. When the star formation activity is poor, the value of the Q(rJgi) goes down.

For each studied galaxy and at each wavelength, the spiral arms are first unwrapped (Iye et al. 1982), by plotting them in an ln  R versus θmap (see Figures 1 and 2). All the spiral arms in the disk galaxy are inspected in search of color gradients. If a candidate region is found, the arm is straightened, by adding a different phase shift to θ at each value of ln  R, until the arm appears as a horizontal line (see Figure 3). This procedure allows us to add the data from different ln  R, and thus to increase the signal-to-noise ratio of the light profile as a function of θ. One-dimensional plots of Q versus θ are obtained from the unwrapped and straightened galaxy images. These are then compared to one-dimensional plots of (g − J) that trace the dust lane location, and to one-dimensional plots of the K-band (or K' or K_s, depending on what we have available) data that trace the density wave (under the assumption that the near-IR emission of red supergiants can be neglected). The regions that seem to match theoretical expectations according to these plots are selected and compared to stellar population synthesis models, that provide Q as a function of stellar age, t_age.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

If one assumes that stars form in the site of the shock, and that they age as they move away from this location, then distance from the dust lane (at constant radius) parameterizes stellar age. In fact, stretching the model Q(t_age) to fit the data (where Q is a function of angular or linear distance) fixes the ratio between the distance from the shock, d, and the age of the stellar population. If, in addition, the rotational velocity is known, it is possible to find the angular velocity of the spiral pattern, Ω_p, as follows:

$$\begin{equation} \Omega _{p} \cong \frac{1}{R_{\rm mean}} \left(v_{\rm rot}(R) - \frac{d}{t_{\rm age}} \right). \end{equation} \tag{ 3 }$$

Here, R_mean is the mean radius of the region where the gradient has been detected, measured from the center of the galaxy; v_rot(R) is the rotational velocity of the disk at such radius (obtained from the literature); d is the distance from the shock, and t_age is the stellar age obtained from fitting the population synthesis model to the observations. With the pattern speed, we obtain the location of major resonances: the outer Lindblad resonance (OLR), corotation, and the 4/1 resonance. This resonance positions are compared with the observed spiral end points.

As mentioned before, this procedure implicitly assumes circular motion for the stars involved, even though no stellar orbits have been computed.

Shock solutions obtained from semianalytical approaches have shown that gas streamlines are not circular. Roberts (1969), for example, found that streamlines appear as sharp-pointed ovals (see Figure 4). Although such solutions were found for the tightly wound approximation, they provide a good reference for the study of real galaxies (see also Saaf 1974). Here, we obtain gas streamlines under the assumption that young stars follow the motion of the molecular gas, and by adopting a procedure similar to that of Gittins & Clarke (2004).⁴ Our adopted model has a pattern speed of 13 km s⁻¹ kpc⁻¹, and a flat rotation curve with v_rot = 220 km s⁻¹. From R = 5 kpc to R = 10 kpc, streamlines were calculated in steps of 0.05 kpc. Since the stellar models give us Q as a function of t_age, we can assign a stellar age to each point in the streamlines; in order to do so, we need to assume that the onset of star formation occurs at a certain orbital time t, where t = 0 corresponds to the shock position.

Zoom In Zoom Out Reset image size

3.1. Star Formation Onset Delay

The gas simulations were performed with a version of the zeus code (Stone & Norman 1992a, 1992b), that is an Eulerian, time-explicit, finite-difference code for ideal MHD simulations. We employ a two-dimensional grid with 500 × 500 points in polar coordinates. The radial extent goes from R∼ 3 to 30 kpc, and the azimuthal one from θ ∼ 0 to π radians. We assume that the half-disk simulation data have a 180° "mirror" symmetry. The gaseous disk follows an isothermal equation of state with a temperature of 10900 K. At the beginning of the simulation, the magnetic field has a toroidal geometry and a value of 5μG at 8.5 kpc from the galactic center, although it rapidly evolves away from this setup. No self-gravity is included.

The background gravitational potential has two components: one axisymmetric, that consists of a bulge, a disk, and a halo, while the other is the non-axisymmetric spiral arm perturbation. The adopted bulge and disk are described in Allen & Santillan (1991). The dark matter halo is NFW type (Navarro et al. 1996, 1997), with density

$$\begin{equation} \rho (R) = \frac{\rho _{h}}{ {R}/{a_{h}} \left(1 + {R}/{a_{h}} \right)^{2}}, \end{equation} \tag{ 4 }$$

where ρ_h = 1.021 × 10⁻² $M_{\odot }\rm \;pc^{-3}$ and a_h = 15.133 kpc. The two-arm spiral potential has a pitch angle of 155, and is self-consistent in the stellar orbits sense (Pichardo et al. 2003). The simulation is performed in this stellar arms reference frame that rotates with velocity Ω_p = 20 km s⁻¹ kpc⁻¹.

We allow the simulation to evolve during 8 × 10⁸ yr (see Figure 5). By this time, the gas develops spiral arms as a response to the imposed perturbation, initial transients subside, and changes in the simulation are only observable at very long timescales. Therefore, we assume that density and velocity distributions at this time represent a close approximation to a steady state solution.

Zoom In Zoom Out Reset image size

As shown in Figure 5, the gas in the simulation has responded to the spiral perturbation with a four-arm pattern (Martos et al. 2004). A hypothetical circular gas orbit would encounter four shocks before completing its transit around the disk. In this case, the corresponding velocity field would be very different from the case where only two shocks are considered. For the present investigation we focus on the main shock (i.e., the one closer to the potential minimum of the spiral perturbation), although spiral-arm triggering of star formation may also take place in the secondary shocks.

In order to obtain the orbit a parcel of gas would follow, we locate the main shock position by searching for the gas density maxima closest to the (stellar) spiral arms. These are taken as the start point for the integration of the gaseous orbits. The orbits are then calculated from the velocity data (in the spiral pattern reference frame), using a Runge–Kutta method, every ΔR = 0.005 kpc, and Δt = 7 × 10¹² s. A sample of orbits is shown in Figure 6. Other gas orbits with MHD simulations using the full time evolution of the gas can be seen in Gómez & Martos (2009), who find a close resemblance between gas orbits and the central family of ballistic orbits in stellar dynamics.

Zoom In Zoom Out Reset image size

5.1. Q Profiles

From the first approach, involving semianalytical solutions, we obtain Q index profiles for the non-circular and circular cases, respectively, calculated in a narrow annulus with mean radius ∼7 kpc; we show them in Figure 7. To produce the Q profiles, we assume a background population of old stars with an age of 5 × 10⁹ yr (i.e., the density wave), plus a young burst of star formation lasting 2 × 10⁷ yr. Both the old and the young populations have a Salpeter IMF, with M_lower = 0.1 M_☉ and M_upper = 10 M_☉.⁶ For the present case, we adopt a constant fraction of young stars of 2% by mass (see Equation (9) in MG09).⁷ We then use the stellar population synthesis models by S. Charlot & G. Bruzual (2007, private communication) to find the changing g, r, i, and J emission of the population mixture as it evolves.

Zoom In Zoom Out Reset image size

In the circular case, the stellar complexes that give origin to the Q profile were all born at the same (ln R, θ) position. In the non-circular case, due to the trajectories of the gas streamlines, the stars that give rise to the Q profiles were not all formed at the same location. The age of the stellar complexes that contribute to the Q profile as a function of distance from the shock position is plotted in Figure 8. As expected, young objects tend to spend more time concentrated toward the spiral shock (d = 0 kpc), with the consequence that the peak of the Q profile (dashed line in Figure 7) occurs closer to the shock.

Zoom In Zoom Out Reset image size

Another aspect to be noticed in Figure 7 is that there is no "downstream fall" of the gradients (i.e., Q does not fall below, in this case, ∼1.56, contrariwise to what is sometimes observed with the data). Martínez-García et al. (2009) had hypothesized that the fall of the observed Q profiles below the model values (assuming pure circular orbits) might be caused by stellar non-circular motions in the data. However, results may differ if ballistic trajectories with post-shock velocities are considered (e.g., Bash et al. 1977; Bash 1979).

5.2. Pattern Speeds

From the Q profiles, it is possible to derive pattern speeds by comparing stellar population synthesis models (Q versus t_age) to the "data" (Q versus d), assuming implicitly that stars move in circular orbits, as was done by GG96 and MG09, and was sketched above in Section 2.

In the case of the semianalytical approach, "data" Q profiles were obtained from narrow annular regions, ∼0.001 kpc wide. Given that we know the input pattern speed, we can search for systematic effects in the determination of Ω_p. Figure 9 shows the input pattern speed (long-dashed line), and Ω_p values that are obtained from the color gradients for three different presuppositions about the delay of star formation onset after the shock (0, 10, and 20 Myr). For reasons explained below in Section 6, rather than stretching the stellar population synthesis model to the "data," we just compare the positions of the maxima of the model Q, on the one hand, and of the shock location, on the other.

Zoom In Zoom Out Reset image size

There is a significant systematic effect, whereby the derived Ω_p is always larger than the input pattern speed; moreover, the effect decreases with galactocentric radius. Also, the difference between the input and the output pattern speeds decreases for larger star formation onset delays.

In the case of the MHD simulations, we again analyze the pattern speeds derived when assuming a circular motion dynamic model. In order to obtain Ω_p, we compare the peak value of Q in the stellar population model to the density maximum of the gas in the simulations. In Figure 10, we show the Ω_p values obtained for δt_shock of 0, 10, and 20 Myr.

Zoom In Zoom Out Reset image size

One of the caveats of the GG96 method, that was not discussed in MG09, comes from image processing. Each pixel in the optical and infrared images actually contains information about many young stellar orbits that all fall within the same spatial region. The process of "unwrapping" the spiral arms then averages pixels in the θ and ln  R directions (the "straightening" of the arms is just a shift in the θ direction and does not involve pixel averaging). In order to increase the signal-to-noise ratio of the Q profiles, an additional averaging in the ln  R direction is done in selected regions. This last step averages together orbits with different angular speeds, and hence stars of slightly different ages at a fixed distance from the shock.

Figure 11 shows the behavior of Q profiles from an unwrapped and straightened synthetic image (see below), under the assumption of pure circular motion; the different profiles are taken from various mean radii, after averaging in ln  R. The shift of the peak toward d = 0 is expected as we approach the corotation radius. However, due to the pixel averaging, a drop in the maximum value of Q is also obtained. If, in order to match the observed or synthetic profile peak, a constant downward shift is applied to the stellar population synthesis model,⁸ the "wings" of the profiles will not match, since those of the stellar model will now lie within the "wings" of the data. This mimics wider gradients, and in order to fit the observed or synthetic profile we will need to overstretch the model. Consequently, the pattern speed obtained via Equation (3) will be lower, because the stellar age (t_age) at a fixed distance d from the shock will be underestimated. The effect amounts to ∼1 to 2 km s⁻¹ kpc⁻¹ for moderate pixel averaging, but can be as large as 10 km s⁻¹ kpc⁻¹ when the averaged region in ln  R is very extended.

Zoom In Zoom Out Reset image size

One way out of this problem is to compare only the distances between the onset of star formation and the Q peak, in the model and the data, instead of the whole Q profiles, since image processing affects less the peak positions than the profile shapes. We have used this approach with our semianalytical calculations and MHD simulations; unfortunately, in the case of real data, the position of the star formation onset is generally unknown (see also Section 3.1).

To better estimate the effects of averaging, with the R, Θ, and t (time elapsed since the shock, t = 0 corresponds to the shock location) data of the streamlines obtained from the semianalytical solutions, we constructed images of a synthetic "galaxy" as would be observed in the g, r, i, and J bands with the aid of stellar population synthesis models (S. Charlot & G. Bruzual 2007, private communication; the fraction of young stars constitute is 2% by mass). The model face-on galaxy is at a distance of 35 Mpc and the synthetic images have 1'' × 1'' pixels. Whenever various streamlines were located at the same pixel, their elapsed times since the shock were averaged. Next, the synthetic images were unwrapped. Each hyperpixel of the unwrapped images is 3° long in the θ-axis; the ln  R-axis has 250 hyperpixels in total, from the center of the galaxy to R = 10 kpc. Then, the synthetic unwrapped images were straightened according to the shock's pitch angle (∼5266) that is smaller than the potential minimum's pitch angle (∼5739; see the Appendix). Finally, Q profiles were obtained from these "images" in sections with a width of five hypercolumns in ln  R. Pattern speeds were measured by comparing the "observed" Q (versus d) from the simulated images, at different radii, with Q (versus t_age) from stellar population synthesis models.⁹

In Figure 12, we illustrate the effects of averaging. The dotted line is the input pattern speed. If the synthetic galaxy has circular orbits and the galaxy is analyzed assuming circular orbits, then one obtains the Ω_p values depicted by squares. The solid black squares are obtained with the GG96 "stretching" method, whereas the comparison between model and data Q maxima positions relative to the shock yields the open symbols. Clearly, the solid black squares are biased while the open ones are not.

Zoom In Zoom Out Reset image size

On the other hand, if the synthetic galaxy has non-circular orbits but is analyzed under the assumption of circular orbits, then one obtains the Ω_p values shown as triangles. The solid black triangles are the values obtained with the GG96 method; the empty triangles are the pattern speeds found when comparing only the positions, relative to the shock location, of model and "observed" Q peaks; and the solid line represents the biased values measured from the semianalytical solutions, without averaging, if existing non-circular motions are ignored. The reason for this bias toward higher measured pattern speeds is that real gas streamlines in a steady rotating spiral shock turn somewhat along the arm after passing through the shock. Consequently, stars close to the shock are slightly older than would be expected in a circular model with the same Ω_p. The observer using a circular model would infer that the gas flow into the arms is smaller than is actually the case, i.e., that the difference between the stellar orbital velocity and the pattern speed is smaller than in reality; inward of corotation, this means that the observer would overestimate the pattern speed. Also inside corotation, the effect decreases with galactocentric radius, as the shock strength diminishes.

The open triangles product of the peak comparison are much less sensitive to the competing systematics introduced by orbit averaging. However, ironically, the original "stretching" method (solid black triangles) seems much less biased toward higher values if non-circular motions are neglected, although the dependence of the bias on galactocentric radius is still noticeable. This is a fortunate turn of events, because (1) in real images it is hard to pinpoint the location of the star formation onset and (2) the radial dependence of the bias allows us to detect the presence of non-circular motions and to confirm the link between star formation and disk dynamics.

From theory (e.g., Gittins & Clarke 2004), it is well known that the spiral shock strength diminishes as one approaches corotation. The shock weakening implies fewer radial movements and more circular trajectories. Figures 9 (corotation radius, R_CR, located at ∼17 kpc) and 10 (R_CR ∼ 11 kpc) support this conclusion; they also allow us to predict that, in general, Ω_p values derived from color gradients observed inside corotation in real data will be too high if a circular dynamic model is adopted.¹⁰ For color gradients observed beyond corotation, we expect again to overestimate the pattern speed, but the difference between the actual and the measured Ω_p will now grow with radius. In this case, the radial velocities of the streamlines after the shock have an outward, rather than an inward, component (see the leading case in Roberts 1969). If, however, spiral shocks are not as strong as inside corotation, inverse post-corotation gradients will yield unbiased pattern speed measurements.

How would these results affect the Martínez-García et al. (2009) conclusions about the end points of spiral patterns? Their Figure 4 displays the ratio, $R_{\rm {res}}{/}R_{\rm {end}}$, of the resonance radii (calculated from their pattern speed measurements and rotation curves) to the observed spiral end points in their deprojected, near-infrared galaxy images. If their Ω_p measurements (most of them inside corotation) are systematically too high, then the resonance radii and their ratio to the spiral end points would be systematically smaller than the true values.

In order to better quantify the impact of non-circular motions on the MG09 results, we define the expression $\delta \Omega _{p} = (\Omega _{\rm {data}}- \Omega _{p}^{\prime }) / \Omega _{p}^{\prime }$. Here, $\Omega _{\rm {data}}$ is the pattern speed obtained from the data, by comparing the observed color gradient candidates with the stellar models, under the assumption that stars move in purely circular orbits; Ω'_p is the pattern speed that the spiral should have if it ends at certain resonance. For the 4:1 resonance, for example, we have

$$\begin{equation} \Omega _{p}^{\prime } \sim \frac{v_{\rm {rot}}}{R_{4:1}} \left(1-\frac{\sqrt{2}}{4} \right), \end{equation} \tag{ 5 }$$

where $v_{\rm {rot}}$ is the circular orbital velocity (obtained from the literature, see MG09). And if the spiral ends at the OLR, then

$$\begin{equation} \Omega _{p}^{\prime } \sim \frac{v_{\rm {rot}}}{R_{\rm OLR}} \left(1+\frac{\sqrt{2}}{2} \right). \end{equation} \tag{ 6 }$$

To calculate δΩ_p for the 23 regions analyzed in MG09, we take both $\Omega _{\rm {data}}$ and $R_{\rm {end}}$ from their Table 4. In Figures 13 and 14, we plot δΩ_p versus $R_{\rm {mean}}{/}R_{\rm {end}}$, where $R_{\rm {mean}}$ is the mean orbital radius of the studied region. Figure 13 corresponds to the case where the end point of the spiral pattern is fixed at the 4:1 resonance. Contopoulos & Grosbol (1986) state that strong spirals (Hubble types Sb or Sc) may be truncated at this resonance. The δΩ_p values show the trend with radius expected from Figures 9 and 10. This would imply that we are actually detecting non-circular motions from the photometric data!

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

However, three discrepancies with respect to theoretical expectations are also present. The first one is related to the systematic effects discussed in Section 6. In spite of the tendency to actually measure values around the true pattern speed with the method used by GG96 and MG09 (solid black triangles in Figure 12), there is a lack of points below the δΩ_p = 0 line. The second difference is the gap around $R_{\rm {mean}}/R_{\rm {end}} \sim 0.6$, where a point with huge vertical error bars is located.¹¹ If spirals are really truncated at the 4:1 resonance, and hence inside or near corotation, one would expect a continuous distribution of points between the inner Lindblad resonance and the spiral end point (i.e., $R_{\rm {mean}}/R_{\rm {end}}=1$). The third inconsistency concerns the actual values attained by δΩ_p at low $R_{\rm {mean}}/R_{\rm {end}}$. From the solid black triangles in Figure 12, the maximum expected value of δΩ_p is ≲0.5. The high δΩ_p values would imply the presence of extremely strong shocks near the spiral arms.

Figure 14, in contrast, shows the behavior of δΩ_p for the case when the spiral end points are fixed at the OLR. If spirals stop near this resonance, the expected $R_{\rm {mean}}{/}R_{\rm {end}}$ value for corotation is ∼0.59, and the gap in the data distribution is naturally explained. We also observe now that some points are located below the δΩ_p = 0 line, as expected from Figure 12, owing to the systematic error introduced by the "stretching" GG96 method and discussed in Section 6. The three last points, marked with a long-dashed ellipse, represent the "inverse color gradients" (i.e., located beyond corotation) examined in MG09 (these points were treated as being inside corotation in Figure 13). The points surrounded by a short-dashed ellipse belong to the galaxy NGC 578. If NGC 578 actually ends before or at corotation (as argued by MG09), its position in Figure 13 is in concordance with theoretical expectations (i.e., close to the δΩ_p = 0 line). On the other hand, NGC 578 may be the only galaxy in the MG09 sample with a spiral pattern that stops at the 4:1 resonance, as predicted by Contopoulos & Grosbol (1986).

Under the assumption that the orbits of young stars preserve the velocity components of the parent molecular clouds where they form, we have analyzed the effects that non-circular motions would have on azimuthal color gradients. Semianalytical calculations and MHD simulations show that the spiral pattern speeds derived from the comparison between color gradients and stellar population synthesis models, assuming purely circular motions, would have values systematically higher than the real ones for regions within corotation. Also inside corotation, the effect would decrease with galactocentric radius. Non-circular motions, however, would not prevent the detection of legitimate azimuthal color gradients in real galaxies.

Using a synthetic image, we have also analyzed the effects of image processing and pixel averaging on the method applied in GG06 and MG09 to detect color gradients and derive pattern speeds. We have found that pixel averaging (due to image processing) systematically decreases the derived Ω_p, such that it nearly compensates for the systematic effect introduced when neglecting non-circular motions in the analysis. The net result is that the correct spiral pattern speeds and resonance locations can be obtained. Nevertheless, a residual trend of slightly higher pattern speeds at lower radii can still be discerned (solid triangles in Figures 12), so that it is possible to detect the presence of non-circular motions and confirm the link between star formation and disk dynamics.

We have re-examined the results obtained by MG09. When normalizing the mean radii where gradients were found by the end point radii, in order to treat the whole sample as a single galaxy, we were able to reproduce the trend of δΩ_p with radius expected if non-circular motions are ignored (as these authors did). The size of the observed δΩ_p in the data only matches the theoretical expectations if spiral patterns end at the OLR, and not at the 4:1 resonance (compare Figure 12 with both Figures 13 and 14).

Future studies of azimuthal color gradients in disk galaxies by orbit calculations should include the effects of disk heating (Asiain et al. 1999) and a correct modeling of the IMF as a discrete distribution. These factors may be important, since high-mass stars may take different trajectories when compared to low-mass stars, due to the effects of dynamical friction (Chandrasekhar 1943).

We thank the anonymous referee for his/her very constructive suggestions and comments that helped us improve the presentation of this work.

E. Martínez-García acknowledges financial support from CIDA (Centro de Investigaciones de Astronomía, Venezuela). Special thanks to César Briceño, Katherina Vivas, Carlos Abad, Gustavo Bruzual, Gladis Magris, Jesús Hernández, and Francisco Fuenmayor, also at CIDA. G.C.G. thanks DGAPA (UNAM) grant IN-106809 and CONACyT grant 50402-F. R.A.G.L. acknowledges DGAPA (UNAM) grant IN111007 and CONACyT grant 48589-F.

The differential equations needed to obtain the density and velocity components, product of spiral shocks, are fully explained in Roberts (1969), Shu et al. (1972, 1973), and Gittins & Clarke (2004).¹² The independent parameters used in this investigation for the solutions are: number of arms, m = 2; spiral arms pitch angle i such that sin i = 0.1; angular speed of the spiral pattern, Ω_p = 13 km s⁻¹ kpc⁻¹; effective speed of sound, a = 8 km s⁻¹; amplitude ratio of the perturbed spiral field to the axisymmetric field, F = 0.05; rotation velocity of material,¹³v_rot = 220 km s⁻¹. In Figure 15, we show the positions of the sonic point and the spiral shock relative to the potential minimum (η = 0), with the parameters described above. As already noticed by Tosa (1973) and Gittins & Clarke (2004), the shock front moves away from the potential minimum as the radius increases. In this case, for R ∼ 5 kpc, the shock front is slightly ahead of the potential minimum. These relative positions may change for another set of parameters.

Zoom In Zoom Out Reset image size

The solutions are obtained in curvilinear coordinates η (perpendicular to the spiral equipotential curves) and ξ (parallel to the spiral equipotential curves). We adopt the η and ξ definitions from Shu et al. (1973) and Gittins & Clarke (2004). These definitions differ from those used in Roberts (1969) and Shu et al. (1972) by a multiplicative factor m/sin i, and by an additive constant, chosen such that η = 0 determines the location of the potential minimum. The corresponding expressions in terms of the radius, R, and the angular coordinate in the rotating frame of the spiral pattern, Θ = θ − Ω_pt, are as follows:

$$\begin{equation} \eta = \frac{m}{\tan {i}} \ln \left({\frac{R}{R_{a}}}\right) + m (\theta - \Omega _{p}t) + \pi, \end{equation} \tag{ A1 }$$

$$\begin{equation} \xi = -m \ln \left({\frac{R}{R_{a}}}\right) + \frac{m}{\tan {i}} (\theta - \Omega _{p}t). \end{equation} \tag{ A2 }$$

The η coordinate is obtained from the shock solution, together with the velocity components u_ξ = u_ξ0 + u_ξ1 and u_η = u_η0 + u_η1. The subscript 0 labels the unperturbed velocity (i.e., in the absence of a spiral perturbation), and the subscript 1, the perturbation due to the spiral gravitational field. The ξ coordinate is obtained by solving the equation

$$\begin{equation} \frac{d\xi }{d\eta }=\frac{u_{\xi }}{u_{\eta }} \end{equation} \tag{ A3 }$$

that can also be expressed as

$$\begin{equation} \int ^{\eta _{\rm stream}}_{\eta _{\rm sub}} \left(\frac{d\xi }{d\eta }\right) d\eta = \xi _{\rm stream} - \xi _{\rm sub}. \end{equation} \tag{ A4 }$$

The shock on the subsonic branch occurs at η_sub and ξ_sub. η_stream and ξ_stream are coordinate values at any point along the streamline; their maximum possible values are η_sup and ξ_sup, i.e., the point at which the shock occurs on the supersonic branch.¹⁴

In order to obtain ξ_sub, we proceed in the way described below. We define the quantity

$$\begin{equation} x = R_{0} \left({ \frac{\int ^{\Theta _{\rm sup}^{\prime }}_{\Theta _{\rm sub}^{\prime }} \left(\frac{R}{R_{a}}\right)^{\prime } d\Theta ^{\prime }}{\Theta _{\rm sup}^{\prime } - \Theta _{\rm sub}^{\prime }} } \right)^{-1}, \end{equation} \tag{ A5 }$$

where $(\frac{R}{R_{a}})^{\prime }$ and Θ' are the $(\frac{R}{R_{a}})$ and Θ values obtained when ξ_sub = 0 and R_a = 1. R₀ is the average radius of the streamline.

We use the x value to get

$$\begin{equation} \xi _{\rm sub} = {-}\!\left(m+\frac{m}{\tan ^{2}{i}}\right) \ln {x}. \end{equation} \tag{ A6 }$$

With this definition, the mean value of $(\frac{R}{R_{a}})$ will be R₀.

The time elapsed since the shock, t, is obtained by integrating the velocity along the Θ direction

$$\begin{equation} t = \int ^{\Theta _{\rm stream}}_{\Theta _{\rm sub}} \frac{R d\Theta }{u_{\Theta }}, \end{equation} \tag{ A7 }$$

where t = 0 when Θ_stream = Θ_sub, and

$$\begin{equation} u_{\Theta } = \frac{u_{\eta } + u_{\xi }\cot {i}}{\cot {i}\cos {i} + \sin {i}}. \end{equation} \tag{ A8 }$$