The Kugel-Khomskii (KK) Hamiltonian describes spin and orbital superexchange interactions between d1 ions in an ideal cubic perovskite structure, in which the three t2g orbitals are degenerate in energy and electron hopping is constrained by cubic site symmetry. In this paper we implement a variational approach to mean-field theory in which each site i has its own n×n single-site density matrix ρ(i), where n, the number of allowed single-particle states, is 6 (3 orbital times 2 spin states). The variational free energy from this 35 parameter density matrix is shown to exhibit the unusual symmetries noted previously, which lead to a wave-vector-dependent susceptibility for spins in α orbitals which is dispersionless in the qα direction. Thus, for the cubic KK model itself, mean-field theory does not provide wavevector “selection,” in agreement with rigorous symmetry arguments. We consider the effect of including various perturbations. When spin-orbit interactions are introduced, the susceptibility has dispersion in all directions in q space, but the resulting antiferromagnetic mean-field state is degenerate with respect to global rotation of the staggered spin, implying that the spin-wave spectrum is gapless. This possibly surprising conclusion is also consistent with rigorous symmetry arguments. When next-nearest-neighbor hopping is included, staggered moments of all orbitals appear, but the sum of these moments is zero, yielding an exotic state with long-range order without long-range spin order. The effect of a Hund’s rule coupling of sufficient strength is to produce a state with orbital order.