  A mixture of two kinds of identical bosons, species $1$ and species $2$, held
in a harmonic potential and interacting by harmonic intra-species and
inter-species particle-particle interactions is discussed. To prove
Bose-Einstein condensation of the mixture three steps are needed. First, we
integrate the all-particle density matrix, employing a four-parameter matrix
recurrence relations, down to the lowest-order intra-species and inter-species
reduced density matrices of the mixture. Second, the coupled Gross-Pitaevskii
(mean-field) equations of the mixture are solved analytically. Third, we
analyze the mixture's reduced density matrices in the limit of an infinite
number of particles of both species $1$ and $2$ (when the interaction
parameters, i.e., the products of the number of particles times the
intra-species and inter-species interaction strengths, are held fixed) and
prove that: (i) Both species $1$ and $2$ are 100\% condensed; (ii) The
inter-species reduced density matrix per particle is separable and given by the
product of the intra-species reduced density matrices per particle; and (iii)
The mixture's energy per particle, and reduced density matrices and densities
per particle all coincide with the Gross-Pitaevskii quantities. Finally, when
the infinite-particle limit is taken with respect to, say, species $1$ only
(with interaction parameters held fixed) we prove that: (iv) Only species $1$
is 100\% condensed and its reduced density matrix and density per particle, as
well as the mixture's energy per particle, coincide with the Gross-Pitaevskii
quantities of species $1$ alone; and (v) The inter-species reduced density
matrix per particle is nonetheless separable and given by the product of the
intra-species reduced density matrices per particle.
