Q$_{min}$ ∼ 1 . 3 at M$_{d}$ ∼ 4 . 4 × 10 $^{9}$M ⊙ in the lower right. Not surprisingly, M$_{d}$ for scenario KL is much larger, since there is no hole in the mass distribution. Since the radial dispersion σ$_{R}$ and epicycle frequency κ are constrained by the data, the increased surface density yields lower values for Q$_{min}$ . For scenario KL, the best fit region extends below the Q$_{min}$ = 1 line, which suggests that some of the best fit models are locally unstable. The PDF in Q$_{min}$ and X$_{min}$ is found in Fig. 9. For scenario K, this plot shows that X$_{min}$ lies between 1 and 3, and appears to correlate weakly with Q$_{min}$ , in the sense that models with high Q$_{min}$ also have high X$_{min}$ . Fig. 9 also shows that all models from scenario KL have X$_{min}$ < 1, suggesting very strong instability to global modes. To discern whether or not the models for scenario KL are realistic, we must determine whether or not galaxies with Q < 1 are possible. Unstable galaxies generate spiral structure that heat stars as they propagate outwards but continual cooling could occur if gas infall onto the surface of the galaxy were substantial enough. Such gas would contribute to star formation. Fuchs (1999) conducted a maximum disc analysis for NGC 6503 suggesting that Q < 1 implies a star formation rate of 40 M ⊙ yr − $^{1}$; this contradicts the rate of 1.5 M ⊙ yr − 1 obtained from the observed H α flux (Kennicutt et al. 1994). Therefore, it is unlikely that this cooling method is important in NGC 6503. In the absence of star formation, a disc in which Q < 1 may fragment as small scale perturbations grow due to the Toomre instability. The likely evolution for such models is to eventually recombine into a single stable disc, but not before the stellar dispersions have increased dramatically. We select several models from the PDF for both runs (identified by stars on Fig. 8) and evolve them forward in time to examine their stability. To do this, we employ Dehnen’s (2000) N -body algorithm implemented by Stiff (2003), a fast multipole method that produces nearly O ( N ) scaling. All models are evolved for 5 Gyr with a time-step of 0 . 5 Myr and a softening parameter of 25 pc. We use 500K particles for the disc, 50K for the bulge and 1M for the halo. The bar evolution is quantified using the magnitude of the second Fourier mode A , which measures the strength of two-armed asymmetries (e.g., Shen & Sellwood 2004): where θ$_{j}$ is the azimuthal coordinate of the j th particle in the disc. The bar evolution is found in Fig. 10. The results show the regions of the M$_{d}$ − Q$_{min}$ plane that are most susceptible to bar formation. For scenario K, we find that much of the region is at least mildly bar unstable, with the lower right region more strongly unstable (note that the red curves in Fig. 10 correspond to the discs with lowest Q$_{min}$ in each run) and the upper left only showing a mild 16 A = ∣ N ∑ j =1 exp(2 iθ$_{j}$ ) N ∣ (18)