  We construct several new families of exactly and quasi-exactly solvable
BC_N-type Calogero-Sutherland models with internal degrees of freedom. Our
approach is based on the introduction of two new families of Dunkl operators of
B_N type which, together with the original B_N-type Dunkl operators, are shown
to preserve certain polynomial subspaces of finite dimension. We prove that a
wide class of quadratic combinations involving these three sets of Dunkl
operators always yields a spin Calogero-Sutherland model, which is
(quasi-)exactly solvable by construction. We show that all the spin
Calogero-Sutherland models obtainable within this framework can be expressed in
a unified way in terms of a Weierstrass P function with suitable half-periods.
This provides a natural spin counterpart of the well-known general formula for
a scalar completely integrable potential of BC_N type due to Olshanetsky and
Perelomov. As an illustration of our method, we exactly compute several energy
levels and their corresponding wavefunctions of an elliptic quasi-exactly
solvable potential for two and three particles of spin 1/2.
