  We introduce a variant of the multi-grid Monte Carlo (MGMC) method, based on
the embedding of an $XY$ model into the target model, and we study its
mathematical properties for a variety of nonlinear $\sigma$-models. We then
apply the method to the two-dimensional $O(N)$-symmetric nonlinear
$\sigma$-models (also called $N$-vector models) with $N=3,4,8$ and study its
dynamic critical behavior. Using lattices up to $256 \times 256$, we find
dynamic critical exponents $z_{int,{\cal M}^2} \approx 0.70 \pm 0.08$, $0.60
\pm 0.07$, $0.52 \pm 0.10$ for $N=3,4,8$, respectively (subjective 68\%
confidence intervals). Thus, for these asymptotically free models, critical
slowing-down is greatly reduced compared to local algorithms, but not
completely eliminated; and the dynamic critical exponent does apparently vary
with $N$. We also analyze the static data for $N=8$ using a finite-size-scaling
extrapolation method. The correlation length $\xi$ agrees with the four-loop
asymptotic-freedom prediction to within $\approx 1\%$ over the interval $12 <
\xi < 650$.
