Deep defects in wide band gap semiconductors have emerged as leading qubit candidates for realizing quantum sensing and information applications. Due to the spatial localization of the defect states, these deep defects can be considered as artificial atoms/molecules in a solid state matrix. Here we show that unlike single-particle treatments, the multiconfigurational quantum chemistry methods, traditionally reserved for atoms/molecules, accurately describe the many-body characteristics of the electronic states of these defect centers and correctly predict properties that single-particle treatments fail to obtain. We choose the negatively charged nitrogen-vacancy (NV$^-$) center in diamond as the prototype defect to study with these techniques due to its importance for quantum information applications and because its properties are well-known, which makes it an ideal benchmark system. By properly accounting for electron correlations and including spin-orbit coupling and dipolar spin-spin coupling in the quantum chemistry calculations, for the NV$^-$ center in diamond clusters, we are able to: (i) show the correct splitting of the ground (first-excited) triplet state into two levels (four levels), (ii) calculate zero-field splitting values of the ground and excited triplet states, in good agreement with experiment, and (iii) calculate the energy differences between ground and exited spin-triplet and spin-singlet states, as well as their ordering, which are also found to be in good agreement with recent experimental data. The numerical procedure we have developed is general and it can screen other color centers whose properties are not well known but promising for applications.